Mixed-Monotonicity Reachability

• Mixed-monotonicity reachability is a method that decomposes nonlinear maps using functions that are increasing in one argument and decreasing in the other to enclose reachable sets.
• The framework applies to continuous, discrete, disturbed, and stochastic systems, and has been extended to traffic networks and neural network verification.
• Key tools include Jacobian-bound constructions, SDP formulations to reduce conservatism, and embedding systems that provide computationally tractable interval over-approximations.

Searching arXiv for papers on mixed-monotonicity reachability and closely related formulations. Mixed-monotonicity reachability is a family of reachability methods in which a nonlinear map or vector field is represented by a decomposition function that is monotone increasing in one argument and monotone decreasing in another, so that reachable sets can be enclosed by propagating lower and upper extremal trajectories rather than enumerating all trajectories. In its standard form, a function $f$ is mixed monotone if there exists $g$ such that $g(x,x)=f(x)$, $g$ is increasing in its first argument, and decreasing in its second; the resulting bounds $$g(\underline{x},\overline{x})\le f(x)\le g(\overline{x},\underline{x})$$ for $\underline{x}\le x\le \overline{x}$ are the basic mechanism behind interval over-approximation. This framework has developed from order-theoretic analysis of nonlinear systems into a reachability toolbox for continuous-time and discrete-time dynamics, disturbed systems, traffic networks, polynomial systems, neural networks, implicit neural networks, and stochastic systems (Yang et al., 2018, Meyer et al., 2019, Tahir, 2023, Sivaramakrishnan et al., 2024).

1. Formal definition and theoretical basis

The standard decomposition-function definition treats mixed monotonicity as a property of a map $f:\mathcal{X}\to\mathcal{T}$: there exists $g:\mathcal{X}\times\mathcal{X}\to\mathcal{T}$ such that $g(x,x)=f(x)$, $g$ is monotonically increasing in its first argument, and monotonically decreasing in its second argument. For order intervals $g$0, this yields the enclosure

$g$1

which is the central reachability relation for interval propagation (Yang et al., 2018).

A recurring point in the literature is the distinction between mixed monotonicity of a vector field and mixed monotonicity of the associated flow map. Under the assumptions that $g$2 is locally Lipschitz on an open set containing the domain, the system is forward complete, and the domain is positively invariant, mixed monotonicity of the vector field implies mixed monotonicity of the flow map. This makes the framework directly relevant to sampled-data abstraction and time-discretized reachability, because one can transfer a decomposition of a continuous-time vector field to the time-$g$3 flow map (Yang et al., 2018).

Two sufficient-condition lines have been especially influential. The first uses bounded partial derivatives: if the partial derivatives of a differentiable map are bounded on a region, then the map is mixed monotone on that region, with an explicit decomposition that adds linear correction terms to compensate for sign-unstable derivatives. A corollary is that every continuously differentiable map on a compact set is mixed monotone. The second uses bounded variation: scalar functions of bounded variation admit a Jordan decomposition into increasing and decreasing parts, hence a mixed-monotone decomposition $g$4. The literature is explicit, however, that broad sufficient conditions do not by themselves guarantee computationally useful enclosures; the quality of the decomposition remains decisive (Yang et al., 2018).

2. Embedding systems and interval reachability

The operational core of mixed-monotonicity reachability is the embedding system. For a decomposition $g$5, one forms

$g$6

which is monotone with respect to the order $g$7. This higher-dimensional monotone system propagates lower and upper bounds simultaneously, and its trajectories can be projected back to the original state space as interval over-approximations (Coogan et al., 2015, Yang et al., 2018).

For disturbed continuous-time systems $g$8 with $g$9, a deterministic worst-case embedding is obtained by fixing the disturbance endpoints: $g(x,x)=f(x)$0 If the embedding trajectory remains in the state space, then

$g(x,x)=f(x)$1

This is the basic finite-time reachable-set over-approximation result for disturbed mixed-monotone systems (Abate et al., 2020).

The same logic underlies the Matlab toolbox TIRA, which implements continuous-time mixed-monotonicity, sampled-data mixed-monotonicity, and discrete-time mixed-monotonicity. TIRA presents mixed monotonicity as a strict generalization of monotonicity: when the original system is monotone, the mixed-monotone construction recovers the classical tight interval reachability result; when the system is not monotone, the duplicated-system construction remains applicable and returns a conservative interval enclosure (Meyer et al., 2019).

A complementary sampled-data line reinterprets the time-$g(x,x)=f(x)$2 flow map $g(x,x)=f(x)$3 as a discrete-time map and shows that the sign-stable sensitivity method is equivalent to a mixed-monotonicity construction on that map. When sensitivities are merely bounded rather than sign-stable, one still obtains a guaranteed interval bound by evaluating the flow at two corner points and adding compensation terms; the resulting over-approximation scales at worst linearly with the state dimension (Meyer et al., 2018).

3. Constructing decomposition functions and reducing conservatism

In generic nonlinear settings, decomposition functions are often built from bounds on Jacobians or sensitivities. The derivative-bounded construction selects, for each coordinate and each input component, either the lower or upper endpoint according to the sign of a center Jacobian, then adds linear correction coefficients when the derivative interval crosses zero. This mechanism is the basis of several discrete-time, sampled-data, and neural-network reachability methods (Yang et al., 2018, Meyer et al., 2018).

A more algebraic construction is available for univariate polynomial functions. Every univariate polynomial is globally mixed monotone with a polynomial decomposition function obtained from a Gram matrix representation of the derivative. If $g(x,x)=f(x)$4 with $g(x,x)=f(x)$5, then

$g(x,x)=f(x)$6

and $g(x,x)=f(x)$7 is a global mixed-monotone decomposition function for the polynomial. This establishes that mixed monotonicity is not only a property that can be certified locally from Jacobian bounds, but also one that every univariate polynomial enjoys globally, with an explicit polynomial decomposition function that can be computed from a Gram-matrix representation (Tahir, 2023).

Tightness is a central issue because decomposition functions are not unique. In the polynomial setting, the nonuniqueness comes from both the PSD split $g(x,x)=f(x)$8 and the Gram representation $g(x,x)=f(x)$9. The paper formulates an SDP

$g$0

and suggests $g$1. The justification is that the interval width $g$2 is bounded by terms proportional to $g$3 and $g$4, although the paper also notes that the tightest objective may depend on the polynomial (Tahir, 2023).

The comparison with local Jacobian-bound decompositions is one of the strongest practical claims in this line of work. For $g$5, the global polynomial decomposition

$g$6

is reported to be much tighter than the local decomposition $g$7 on $g$8, even though $g$9 is global and evaluable (Tahir, 2023).

Conservatism can also be reduced geometrically rather than algebraically. A linear state transformation $$g(\underline{x},\overline{x})\le f(x)\le g(\overline{x},\underline{x})$$0 yields a transformed system $$g(\underline{x},\overline{x})\le f(x)\le g(\overline{x},\underline{x})$$1. If the transformed system is mixed monotone, reachability can be computed in transformed coordinates and mapped back as a parallelotope $$g(\underline{x},\overline{x})\le f(x)\le g(\overline{x},\underline{x})$$2 in the original space. The resulting enclosures can be tighter when the true reachable geometry is tilted or skewed relative to the coordinate axes, and multiple decompositions or multiple transformations can be combined by intersection (Abate et al., 2020).

4. Traffic flow, disturbed systems, and stochastic extensions

Mixed monotonicity became important in traffic-flow reachability because realistic diverging-junction models with FIFO spillback are generally not monotone under any orthant order. For partial FIFO traffic networks, each turning movement is decomposed as

$$g(\underline{x},\overline{x})\le f(x)\le g(\overline{x},\underline{x})$$3

with a FIFO-susceptible component and a non-FIFO component. The network dynamics

$$g(\underline{x},\overline{x})\le f(x)\le g(\overline{x},\underline{x})$$4

are shown to be mixed monotone under Assumption 1, using an explicit decomposition function $$g(\underline{x},\overline{x})\le f(x)\le g(\overline{x},\underline{x})$$5 that evaluates the FIFO part at a modified state $$g(\underline{x},\overline{x})\le f(x)\le g(\overline{x},\underline{x})$$6 and keeps the non-FIFO part at $$g(\underline{x},\overline{x})\le f(x)\le g(\overline{x},\underline{x})$$7. This places partial FIFO traffic networks inside the mixed-monotone class even though they are generally not monotone in the original state (Coogan et al., 2015).

The same embedding machinery supports invariance analysis. For continuous-time mixed-monotone systems with disturbances, certain equilibria of the deterministic embedding yield robustly forward invariant boxes, and asymptotically stable embedding equilibria yield attractive invariant boxes. The backward-time embedding gives complementary results: it over-approximates finite-time backward reachable sets and can certify that the complement of an inner box is robustly forward invariant (Abate et al., 2020).

A stochastic extension replaces deterministic order by stochastic order. For the discrete-time nonlinear stochastic system

$$g(\underline{x},\overline{x})\le f(x)\le g(\overline{x},\underline{x})$$8

a stochastic mixed-monotone system is defined by distributional consistency and stochastic monotonicity/antitonicity conditions for a decomposition $$g(\underline{x},\overline{x})\le f(x)\le g(\overline{x},\underline{x})$$9. The interval update retains the deterministic mixed-monotone form,

$\underline{x}\le x\le \overline{x}$0

but is paired with a probability recursion

$\underline{x}\le x\le \overline{x}$1

The guarantee is that each interval over-approximates the stochastic reachable set and satisfies $\underline{x}\le x\le \overline{x}$2, although the paper is explicit that $\underline{x}\le x\le \overline{x}$3 may deteriorate over time and may eventually become trivial (Sivaramakrishnan et al., 2024).

5. Neural-network reachability and robustness certification

For feedforward neural networks with uncertain inputs, mixed monotonicity has been adapted as a layerwise reachability method that uses local bounds on activation derivatives, propagates Jacobian interval bounds, applies the mixed-monotonicity bound to every partial network, and intersects the resulting intervals. The method is sound for any Lipschitz-continuous activation function provided the function, its derivative, and the global argmin and argmax of the derivative are known. Because it evaluates all partial-network decompositions, the complexity is $\underline{x}\le x\le \overline{x}$4, and the returned interval is the tightest interval approximation obtainable by that mixed-monotonicity framework. On VNN benchmarks, the method is reported to be typically $\underline{x}\le x\le \overline{x}$5 to $\underline{x}\le x\le \overline{x}$6 times tighter than IBP, while remaining applicable to ReLU, tanh, ELU, and SiLU (Meyer, 2021).

The framework has also been extended to uncertain-parameter networks, in which every layer’s weight matrix and bias vector vary within known bounds. The mixed-monotonicity method treats the input, weights, and biases of each partial network as a joint uncertainty vector, propagates Jacobian bounds with respect to that vector, and intersects all partial-network enclosures. In the reported experiments on randomly generated ReLU networks, the method is tighter than ESIP from $\underline{x}\le x\le \overline{x}$7, and for width $\underline{x}\le x\le \overline{x}$8 and depth $\underline{x}\le x\le \overline{x}$9 the example values are mixed monotonicity width $f:\mathcal{X}\to\mathcal{T}$0, ESIP width $f:\mathcal{X}\to\mathcal{T}$1, mixed monotonicity time $f:\mathcal{X}\to\mathcal{T}$2 s, and ESIP time $f:\mathcal{X}\to\mathcal{T}$3 s; ESIP becomes infeasible at larger depths, whereas mixed monotonicity still runs up to $f:\mathcal{X}\to\mathcal{T}$4 and beyond (Meyer, 2023).

For implicit neural networks, interval reachability is framed through inclusion functions rather than only Lipschitz constants. The embedded INN couples lower and upper hidden states through a mixed-monotone map built from the Metzler and non-Metzler parts of $f:\mathcal{X}\to\mathcal{T}$5, and contraction theory ensures convergence of the embedded fixed-point iteration. The resulting output interval is an inclusion function for the INN input-output map. This mixed-monotone contracting framework recovers IBP as the FFNN special case and is proved to be tighter or equal to the IBP reachable set on weight-tied infinite-layer networks when $f:\mathcal{X}\to\mathcal{T}$6. In robust training on MNIST at perturbation radius $f:\mathcal{X}\to\mathcal{T}$7, INNs trained with mixed monotonicity achieved average certified robust accuracy of $f:\mathcal{X}\to\mathcal{T}$8, compared with $f:\mathcal{X}\to\mathcal{T}$9 for FFNNs trained with IBP (Davydov et al., 2022).

Mixed monotonicity has also been applied to neural ODEs. The method constructs a continuous-time decomposition from Jacobian bounds and a sign-stability shift $g:\mathcal{X}\times\mathcal{X}\to\mathcal{T}$0, then simulates the lifted system once from $g:\mathcal{X}\times\mathcal{X}\to\mathcal{T}$1 to obtain $g:\mathcal{X}\times\mathcal{X}\to\mathcal{T}$2. A distinctive contribution is the use of the homeomorphism property of neural ODE flows to justify a boundary-based method that propagates only the boundary facets of the initial box, unions the resulting enclosures, and then takes their interval hull. The reported comparison is explicitly a trade-off between tightness and efficiency: CORA zonotopes are the tightest, NNV2.0 star sets are generally intermediate, and TIRA interval boxes are the most conservative but the fastest (Sayed et al., 15 Oct 2025).

6. Software, practical trade-offs, and open problems

The software ecosystem reflects the dual character of mixed-monotonicity reachability as both an analytical theory and an implementation pattern. TIRA organizes interval reachability around a hub function and includes contraction/growth bounds, continuous-time mixed-monotonicity, sampled-data mixed-monotonicity, and discrete-time mixed-monotonicity. Its stated emphasis is simplicity and scalability rather than the accuracy of the over-approximations, and the paper explicitly contrasts interval methods with geometric set representations such as zonotopes and polytopes (Meyer et al., 2019).

The JAX toolbox immrax recasts interval analysis and mixed-monotone reachability as function transforms. It implements inclusion-function transforms natif, jacif, and mjacif, and embedding transforms ifemb, natemb, jacemb, and mjacemb. Because intervals are represented as JAX PyTrees, the toolbox inherits Just-In-Time Compilation, GPU acceleration, and Automatic Differentiability. The mathematical core is an inclusion-function-based embedding whose single interval-valued trajectory bounds all trajectories from the original interval initial conditions and disturbances, and the differentiability of interval-derived objectives enables robust control synthesis inside gradient-based optimization loops (Harapanahalli et al., 2024).

Several limitations recur across the literature. Decomposition functions are not unique, and the “best” objective for tightness is not universal; globality does not guarantee optimality, and local Jacobian bounds may still be preferable in some settings. Interval boxes are computationally attractive but geometrically restrictive, so they can be substantially more conservative than zonotopes, star sets, or other richer set representations. In stochastic reachability, the probability bound may deteriorate rapidly with the horizon. In neural-network verification, tightness often improves by intersecting many partial-network decompositions, but this raises computational cost, as in the $g:\mathcal{X}\times\mathcal{X}\to\mathcal{T}$3 partial-network method for FFNNs (Tahir, 2023, Sivaramakrishnan et al., 2024, Meyer, 2021).

The current frontier is therefore not a single theorem but a set of linked technical questions: how to construct tighter decomposition functions, how to choose objectives or coordinate transformations that reduce conservatism, how to extend constructive global results beyond univariate polynomials, and how to preserve tractability when the model class includes disturbances, stochasticity, implicit layers, uncertain parameters, or continuous-depth dynamics. The recent literature suggests that mixed-monotonicity reachability remains attractive precisely where interval geometry, monotone embeddings, and endpoint propagation align with the required scale and certification guarantees, while richer set representations remain preferable when geometric tightness is paramount (Abate et al., 2020, Davydov et al., 2022, Sayed et al., 15 Oct 2025).