Skip to content

Theory

This document proves every mathematical fact the Twindex checker relies on. Each VALID verdict is a consequence of the results below plus exact rational arithmetic; no floating-point computation is involved in checking.

Notation: I = [0, 1], ‖f‖ = max_{x∈I} |f(x)|.

1. The problem

A model fixes a coefficient degree d ≥ 0, a polynomial source g, boundary values α, β (rational; default 0), and a nonlinearity

f(u) = Σ_{m ∈ M} a_m u^m,   every m ∈ M odd,   every a_m ≥ 0.

For θ ∈ ℝ^{d+1} let

q(x; θ) = Σ_{j=0}^{d} θ_j B_j(x),   B_j(x) = C(d, j) x^j (1 − x)^{d−j},

and consider

(P_θ)   −u''(x) + q(x; θ) f(u(x)) = g(x)  on (0, 1),   u(0) = α,  u(1) = β.

The original (v1) model is d = 3, f(u) = u³, g = 1, α = β = 0.

Two facts are used throughout.

  • The B_j are nonnegative on I and sum to 1. Hence θ ≥ 0 (componentwise) implies q ≥ 0 on I. Every prior box Twindex accepts has lo_j ≥ 0.
  • For odd m, define h_m(a, b) = Σ_{i=0}^{m−1} a^i b^{m−1−i}, so that a^m − b^m = (a − b) h_m(a, b). Because t ↦ t^m is nondecreasing for odd m, h_m(a, b) ≥ 0 whenever a ≠ b; by continuity it is also ≥ 0 when a = b. Setting H(a, b) = Σ_m a_m h_m(a, b) ≥ 0 gives f(a) − f(b) = (a − b) H(a, b) with H ≥ 0 continuous.

2. A one-dimensional comparison lemma

Lemma 1. Let c ∈ C(I) with c ≥ 0, and let w ∈ C²(0,1) ∩ C(I) satisfy −w'' + c w ≥ 0 on (0,1) and w(0) ≥ 0, w(1) ≥ 0. Then w ≥ 0 on I.

Proof. Suppose w(x₀) < 0 for some x₀. Let (a, b) be the connected component of the open set {x ∈ I : w(x) < 0} containing x₀. At an endpoint inside (0,1) we have w = 0 by continuity. At an endpoint equal to 0 or 1 we have w ≥ 0 by hypothesis and w ≤ 0 as a limit of negative values, so again w = 0. On (a, b), c w ≤ 0, so −w'' ≥ −c w ≥ 0: w is concave on [a, b]. A concave function vanishing at both ends of an interval is ≥ 0 on it, contradicting w(x₀) < 0. ∎

3. Existence and uniqueness

Theorem 1. If q ∈ C(I) and q ≥ 0, then (P_θ) has exactly one solution u ∈ C²(I).

Proof.

Existence. Let F(t) = Σ_m a_m t^{m+1}/(m+1). Each m + 1 is even, so F is convex and F ≥ 0. Let ℓ(x) = α(1 − x) + βx and write u = ℓ + v with v ∈ H¹₀(0,1). Consider E(v) = ∫₀¹ ( ½ (v' + ℓ')² + q F(ℓ + v) − g (ℓ + v) ) dx. In one dimension H¹₀ ⊂ C(I) with ‖v‖_{L²} ≤ ‖v‖ ≤ ‖v'‖_{L²}, so E is well defined and E(v) ≥ ½‖v'‖²_{L²} − (|ℓ'| + ‖g‖_{L²})‖v'‖_{L²} − C for a constant C, which is coercive. E is strictly convex (strictly convex ½‖v'‖², plus the convex ∫ q F(v) since q ≥ 0, plus a linear term). It is weakly lower semicontinuous: the first term is a norm squared, and the other two are continuous under the compact embedding H¹₀ ⊂ C(I). So E has a unique minimiser u, whose Euler–Lagrange equation is ∫ u' φ' + q f(u) φ − g φ = 0 for all φ ∈ H¹₀, with u = ℓ + v. Since q f(u) − g ∈ C(I), the weak second derivative u'' = q f(u) − g is continuous, and u ∈ C²(I) solves (P_θ) classically, with u(0) = α, u(1) = β.

Uniqueness. If u₁, u₂ are classical solutions, e = u₁ − u₂ satisfies −e'' + c e = 0 with c = q H(u₁, u₂) ≥ 0 (§1), and e(0) = e(1) = 0. Lemma 1 applied to e and to −e gives e = 0. ∎

Special case (v1 model, f(u) = u³, g = 1). Here also 0 ≤ u ≤ x(1 − x)/2. Lower bound: at a negative interior minimum x₀, −u''(x₀) + q u(x₀)³ ≤ 0 < 1, a contradiction. Upper bound: v = x(1−x)/2 − u satisfies −v'' = q u³ ≥ 0 with zero boundary values, so v ≥ 0 by Lemma 1. The checker does not use these bounds.

4. The barrier bound (from a candidate to the true solution)

Theorem 2. Let q ≥ 0 on I and let u solve (P_θ). Let p be a polynomial with p(0) = α and p(1) = β, set r = −p'' + q f(p) − g, and suppose ‖r‖ ≤ ρ. Then for every x ∈ I

|u(x) − p(x)| ≤ ρ · x(1 − x) / 2.

Proof. Put e = u − p and c = q H(u, p). By §1, c ≥ 0, and it is continuous. Subtracting the equations, and using f(u) − f(p) = (u − p) H(u, p),

−e'' + c e = (−u'' + q f(u)) − (−p'' + q f(p)) = g − (r + g) = −r.

Let b = ρ x(1−x)/2, so −b'' = ρ and b ≥ 0. For w± = b ∓ e,

−w±'' + c w± = ρ + c b ± r ≥ ρ − |r| ≥ 0,

and w±(0) = w±(1) = 0, because e = u − p vanishes at both ends. Lemma 1 gives w± ≥ 0, i.e. |e| ≤ b. ∎

Notes. * The bound needs no assumption about how p was obtained: no convergence of a numerical method, no spectral tail estimate and no inverse-operator norm. It also holds pointwise; the checker uses the pointwise value x_i(1−x_i)/2 at each sensor, which is never larger than 1/8. * The hypotheses the checker must verify are exactly: p(0) = 0, p(1) = 0 (exact rational equality), θ ≥ 0 (via the prior box), and ‖r‖ ≤ ρ (§6). With inhomogeneous boundary values the first two become p(0) = α, p(1) = β.

5. Uniform bounds over a parameter box

Theorem 3. Fix p. For θ, c ∈ ℝ^{d+1}, r(x; θ) = r(x; c) + Σ_j (θ_j − c_j) B_j(x) f(p(x)). Consequently, if |θ_j − c_j| ≤ h_j for all j,

‖r(·; θ)‖ ≤ ‖r(·; c)‖ + Σ_j h_j ‖B_j f(p)‖.

Proof. q is linear in θ and r depends on θ only through q f(p). The inequality is the triangle inequality. ∎

So one candidate p and one exact computation give a residual bound ρ_B valid for every θ in a box B = c ± h. If B lies in a prior box with lo ≥ 0, Theorem 2 then applies to every θ ∈ B with the same ρ_B.

Theorem 3′ (affine candidate families). Fix polynomials p₀, s₀, …, s_d with p₀(0) = α, p₀(1) = β and every s_j(0) = s_j(1) = 0, a centre c and half-widths h ≥ 0. For Δ ∈ D = Π_j [−h_j, h_j] put

p_Δ = p₀ + Σ_j Δ_j s_j,     θ = c + Δ.

The residual r(x; Δ) = −p_Δ'' + q(x; c + Δ) f(p_Δ) − g is a polynomial in (x, Δ). Write it as Σ_α Δ^α R_α(x), a finite sum over multi-indices α. Then:

  1. sup_{x∈I, Δ∈D} |r(x; Δ)| ≤ ρ_D := Σ_α h^α ‖R_α‖, where h^α = Π_j h_j^{α_j};
  2. if c + D ⊆ [0, ∞)^{d+1}, then for every θ = c + Δ with Δ ∈ D, |u_θ(x) − p_Δ(x)| ≤ ρ_D x(1−x)/2;
  3. at each x, {p_Δ(x) : Δ ∈ D} = [p₀(x) − Σ_j h_j |s_j(x)|, p₀(x) + Σ_j h_j |s_j(x)|].

Proof. (1) |Δ^α| ≤ h^α on D, then the triangle inequality. (2) For a fixed Δ, p_Δ is a polynomial with p_Δ(0) = α, p_Δ(1) = β, and its residual is bounded by ρ_D by (1). Since θ ≥ 0, q ≥ 0, so Theorem 2 applies to the candidate p_Δ. (3) p_Δ(x) is affine in Δ, so over the box it attains exactly the interval p₀(x) ± Σ_j h_j |s_j(x)|. ∎

With all s_j = 0 this is exactly Theorem 3. The gain comes from choosing s_j close to the sensitivities ∂u/∂θ_j: the first-order terms R_{e_j} then nearly cancel, leaving a residual that is O(|h|²). The checker expands r exactly (rational arithmetic, at most 20 000 monomials in Δ) and computes every ‖R_α‖ with the certified bound of §6.

Sensor tests over the box (§7) use the exact interval from (3). Every θ ∈ c + D is compatible if max(|lo_i − y_i|, |hi_i − y_i|) + s_i ≤ ε at every sensor. No θ is compatible if lo_i − y_i > ε + s_i or y_i − hi_i > ε + s_i at some sensor. Here s_i = ρ_D x_i(1−x_i)/2.

5.1 Monotone two-corner boxes

Lemma 2 (positivity). If q ≥ 0, g ≥ 0 on I and α, β ≥ 0, then the solution of (P_θ) satisfies u ≥ 0.

Proof. Write f(u) = u φ(u) with φ(u) = Σ_m a_m u^{m−1}. Every m − 1 is even and every a_m ≥ 0, so φ ≥ 0, and φ(u(x)) is continuous. Then −u'' + c u = g ≥ 0 with c = q φ(u) ≥ 0, and u(0) = α ≥ 0, u(1) = β ≥ 0. Lemma 1 gives u ≥ 0. ∎

Theorem 6 (monotone dependence on θ). Under the hypotheses of Lemma 2, if 0 ≤ θ ≤ θ' componentwise, then u_θ ≥ u_θ' on I.

Proof. Let w = u_θ − u_θ' and write q = q(·; θ), q' = q(·; θ'). Subtracting the two equations,

−w'' + q (f(u_θ) − f(u_θ')) = (q' − q) f(u_θ').

By §1, f(u_θ) − f(u_θ') = H(u_θ, u_θ') w with H ≥ 0, so the left side is −w'' + c w with c = q H ≥ 0. On the right, q' − q = Σ_j (θ'_j − θ_j) B_j ≥ 0, and f(u_θ') ≥ 0 because u_θ' ≥ 0 (Lemma 2) and f is nondecreasing with f(0) = 0. Since w(0) = w(1) = 0, Lemma 1 gives w ≥ 0. ∎

Corollary 7 (two-corner box). Let B = Π_j [lo_j, hi_j] ⊆ [0, ∞)^{d+1} and assume the hypotheses of Lemma 2. Let p_lo and p_hi be candidates for θ = lo and θ = hi with residual bounds ρ_lo, ρ_hi (Theorem 2). Then for every θ ∈ B and every x ∈ I,

p_hi(x) − ρ_hi x(1−x)/2  ≤  u_hi(x)  ≤  u_θ(x)  ≤  u_lo(x)  ≤  p_lo(x) + ρ_lo x(1−x)/2.

So B ⊆ C_y if every sensor interval lies in [y_i − ε, y_i + ε], and B ∩ C_y = ∅ if some sensor interval misses it entirely.

Proof. lo ≤ θ ≤ hi, so Theorem 6 gives the two middle inequalities. Theorem 2 gives the outer ones. ∎

Corollary 7 needs only two ordinary residual bounds per box, independent of the box width. Its sensor intervals enclose the true range [u_hi(x), u_lo(x)] up to the two corner errors, so they are nearly sharp. The checker certifies g ≥ 0 exactly: every Bernstein coefficient of every piece of a dyadic subdivision is ≥ 0, which implies g ≥ 0 (§6, Fact A). It reports INVALID if some attained value is negative, and INCONCLUSIVE if neither is decided.

6. Certified sup norms of polynomials

Let p have degree at most n, and write it in the Bernstein basis of degree n on I:

p(x) = Σ_{i=0}^{n} b_i C(n, i) x^i (1 − x)^{n−i},
b_i = Σ_{k=0}^{i} [C(i, k) / C(n, k)] a_k     (a_k = monomial coefficients).

Fact A (enclosure). The basis functions are nonnegative on I and sum to 1, so p(x) is a convex combination of the b_i and |p(x)| ≤ max_i |b_i|.

Fact B (endpoints). b_0 = p(0) and b_n = p(1).

Fact C (subdivision). The de Casteljau algorithm at t = 1/2 maps the coefficients b of p on [α, β] to the exact Bernstein coefficients of p on [α, m] and [m, β], where m = (α + β)/2. Each step averages adjacent entries; the implementation keeps integer numerators and multiplies the common denominator by 2ⁿ, so the computation is exact.

Algorithm (sup_abs_bound). Start from the piece I. Keep a best attained value L: the maximum of |b_0| and |b_n| over every piece generated so far, each of which is |p| at a point by Fact B. Repeatedly take the piece with the largest enclosure U = max |b_i|. If U ≤ L, stop: by Fact A every remaining piece satisfies |p| ≤ U ≤ L on it. Otherwise, if the piece has reached the depth limit, record U; if not, split it by Fact C. Return the maximum of L and all recorded U.

Theorem 4. The returned value R satisfies R ≥ ‖p‖.

Proof. The pieces form a partition of I into dyadic intervals; every point lies in some piece that was either recorded (bound U ≤ R), discarded when its enclosure fell to U ≤ L ≤ R, or never split because the loop stopped (its enclosure is at most the popped U, which is ≤ L ≤ R). In each case Fact A bounds |p| on that piece by a quantity ≤ R. ∎

L is also a lower bound on ‖p‖, so R − L measures how tight the result is.

7. From sensor inequalities to the claims

Fix a problem instance (sensors x_i, rational data y_i, tolerance ε, prior box Θ ⊆ [0, ∞)^{d+1}, target index k), and let C_y = {θ ∈ Θ : max_i |u_θ(x_i) − y_i| ≤ ε}.

Suppose Theorem 2 (or Theorem 3 with Theorem 2) gives |u_θ(x_i) − p(x_i)| ≤ s_i. By the triangle inequality:

  • if |p(x_i) − y_i| + s_i ≤ ε for every i, then θ ∈ C_y (certified compatible);
  • if |p(x_i) − y_i| − s_i > ε for some i, then θ ∉ C_y (certified incompatible);
  • otherwise these sufficient conditions decide nothing (inconclusive).

For a box B ⊆ Θ with a uniform bound, the same three cases certify B ⊆ C_y, B ∩ C_y = ∅, or nothing.

Corollary (pair certificate). If θ, θ' ∈ C_y are certified compatible, then D_Q(y) := sup_{φ,φ'∈C_y} |φ_k − φ'_k| ≥ |θ_k − θ'_k|. Moreover every estimator T (any function of the data) satisfies

max( |T(y) − θ_k|, |T(y) − θ'_k| ) ≥ |θ_k − θ'_k| / 2,

because both θ and θ' are consistent with the same observed y at noise level ε. This is the triangle inequality. So no method, PINN or otherwise, can recover θ_k from these data with worst-case error below half the certified gap.

Theorem 5 (paving certificate). Let a binary tree be given whose root box is Θ and in which every internal node splits its box B on an axis j at a point m with lo_j(B) < m < hi_j(B) into the two closed boxes {θ ∈ B : θ_j ≤ m} and {θ ∈ B : θ_j ≥ m}. Suppose every leaf labelled inside is certified ⊆ C_y, and every leaf labelled outside is certified ∩ C_y = ∅, by §5 and the case analysis above. Let In be the inside leaves and Out⁺ the inside and unresolved leaves. Then

max_{B∈In} hi_k(B) − min_{B∈In} lo_k(B)  ≤  D_Q(y)  ≤  max_{B∈Out⁺} hi_k(B) − min_{B∈Out⁺} lo_k(B),

where the left side requires In ≠ ∅. If Out⁺ = ∅ then C_y = ∅.

Proof. Coverage: by induction on the tree, the leaves cover Θ, because each split writes its box as the union of its two children. Upper bound: let θ ∈ C_y ⊆ Θ. Then θ lies in some leaf, and that leaf is not an outside leaf (those are disjoint from C_y), so it belongs to Out⁺ and θ_k ∈ [min_{Out⁺} lo_k, max_{Out⁺} hi_k]. Hence any two points of C_y differ in coordinate k by at most the stated width. If Out⁺ = ∅ no leaf can contain θ, so C_y = ∅. Lower bound: pick an inside leaf A attaining the maximum of hi_k and an inside leaf A' attaining the minimum of lo_k. Then choose θ ∈ A with θ_k = hi_k(A) and θ' ∈ A' with θ'_k = lo_k(A'). Both lie in C_y because inside leaves are subsets of C_y. ∎

The checker verifies the split condition at every node, recomputes every leaf box exactly from the root and the splits, re-certifies every inside and outside leaf, and compares the claimed bounds with the two expressions above.

Theorem 8 (oriented paving). Fix a rational origin o, an invertible rational matrix A and a root box R in coordinates φ, and set θ = o + A φ. Suppose R contains the interval hull of A⁻¹(Θ − o), and take a bisection tree over R whose leaves are φ-boxes P. Write θ(P) for the parallelotope o + A P and H(P) for its exact interval hull. Suppose:

  • every leaf labelled outside either satisfies H(P) ∩ Θ = ∅, or carries a candidate family p_Δ = p₀ + Σ_k Δ_k s_k (with θ = o + A(c + Δ), c the centre of P, Δ ∈ P − c) that the sensor tests of Theorem 3′ prove incompatible;
  • every leaf labelled inside satisfies H(P) ⊆ Θ and is proved compatible in the same way.

Let In be the inside leaves and Out⁺ the inside and unresolved leaves with H(P) ∩ Θ ≠ ∅. Then

max_{In} hi_k(H(P)) − min_{In} lo_k(H(P))  ≤  D_Q(y)  ≤
max_{Out⁺} min(hi_k(H(P)), hi_k(Θ)) − min_{Out⁺} max(lo_k(H(P)), lo_k(Θ)).

Proof. Every θ ∈ Θ has φ = A⁻¹(θ − o) ∈ R, hence lies in θ(P) for some leaf P. A leaf with H(P) ∩ Θ = ∅ contains no point of Θ. For θ ∈ θ(P) ∩ Θ we have θ ≥ 0, so q ≥ 0. Since q(c_θ + AΔ) = q(c_θ) + Σ_k Δ_k E_k with E_k = Σ_j A_{jk} B_j, the expansion of Theorem 3′ with directions E_k bounds the residual of the family member at Δ = φ − c, and Theorem 2 applies. Hence outside leaves contain no point of C_y, and inside leaves (which lie in Θ) contain only points of C_y. The upper bound follows as in Theorem 5, using H(P) ∩ Θ as an enclosure of each remaining leaf. For the lower bound, θ_k over the parallelotope θ(P) attains both ends of H(P)_k at vertices. ∎

Oriented pavings let the leaves follow a thin, slanted compatible set, for example by choosing the columns of A along the singular directions of the sensitivity matrix. Monotone corners (Corollary 7) are not used in oriented pavings, because the componentwise order is not preserved by A.

7.1 Two space dimensions (unit square)

Let Ω = (0, 1)², q(x, y; θ) = Σ_{j,k} θ_{jk} B_j(x) B_k(y) with θ ≥ 0 (so q ≥ 0), f as in §1, g a polynomial, and consider

(P²_θ)   −Δu + q f(u) = g in Ω,   u = 0 on ∂Ω.

Well-posedness. E(v) = ∫_Ω ½|∇v|² + q F(v) − g v is coercive, strictly convex and weakly lower semicontinuous on H¹₀(Ω). Since H¹₀ ⊂ L^r for every r < ∞ in two dimensions, ∫ q F(v) is finite. The unique minimiser is the unique weak solution u. Monotone truncation (Stampacchia) gives u ∈ L^∞, and interior elliptic regularity makes u continuous, in fact smooth, inside Ω.

Theorem 9 (2-D barrier). Let p be a polynomial vanishing on ∂Ω with |−Δp + q f(p) − g| ≤ ρ on Ω̄. Let w ∈ C²(Ω̄) satisfy −Δw ≥ 1 in Ω and w ≥ 0 on ∂Ω. Then |u − p| ≤ ρ w in Ω.

Proof. e = u − p ∈ H¹₀(Ω) is a weak solution of −Δe + c e = −r with c = q H(u, p) ∈ L^∞, c ≥ 0 (§1). The functions v± = ρw ∓ e satisfy −Δv± + c v± ≥ ρ + cρw ± r ≥ 0 weakly, because ρw ≥ 0 on Ω̄ (by the weak maximum principle, −Δw ≥ 0 and w ≥ 0 on ∂Ω give w ≥ 0), and v± ≥ 0 on ∂Ω. The weak maximum principle for −Δ + c with c ≥ 0 (Gilbarg–Trudinger, Theorem 8.1) gives v± ≥ 0 a.e. Since u and p are continuous in Ω, the bound holds pointwise. ∎

The checker uses, pointwise at each sensor, the minimum of the three barriers x(1−x)/2, y(1−y)/2 and (x(1−x) + y(1−y))/4. Each satisfies −Δw = 1 and vanishes or is nonnegative on ∂Ω, so each bound holds separately. The residual bound is certified by the tensor-product Bernstein enclosure (twindex/poly2.py), whose Facts A–C hold coordinatewise. Corner compatibility: a polynomial vanishing on ∂Ω has Δp = 0 at the corners, so its residual there is −g(corner). If g is nonzero at a corner, the solution is not C² up to that corner and ρ ≥ |g(corner)|. This floor can block certification at small ε, though a large enough ε can still close. For small residuals, use sources that vanish at the corners. "p vanishes on ∂Ω" is checked exactly on the coefficient grid: the first row and column vanish, and the row and column sums vanish.

7.2 Time-dependent problems (heat equation)

Let Q = (0, 1) × (0, 1] (time scaled to [0, 1]), q(x; θ) = Σ_j θ_j B_j(x) with θ ≥ 0, f as in §1, g(x, t) and u₀(x) polynomials with u₀(0) = u₀(1) = 0, and consider

(H_θ)   u_t − u_xx + q f(u) = g in Q,   u(0, t) = u(1, t) = 0,   u(x, 0) = u₀(x).

Well-posedness. f is nondecreasing with f(0) = 0 and q ≥ 0, so the comparison principle gives uniqueness and an a-priori L^∞ bound. Standard theory for semilinear parabolic equations with monotone lower-order terms gives a solution that is C^{2,1} in Q and continuous on Q̄. With the compatibility below it is smooth up to the corners to the order used. See Lieberman, Second Order Parabolic Differential Equations, Ch. 2 and 5.

Theorem 10 (parabolic barrier). Let p(x, t) be a polynomial that matches the data exactly: p(0, t) = p(1, t) = 0 and p(x, 0) = u₀(x). Suppose |p_t − p_xx + q f(p) − g| ≤ ρ on Q̄. Then

|u − p| ≤ ρ · min(t, x(1 − x)/2)   in Q.

Proof. e = u − p vanishes on the parabolic boundary. It satisfies e_t − e_xx + c e = −r, where r is the residual of p and c = q H(u, p) ≥ 0 is bounded (§1). Let w be either t or x(1 − x)/2. Both satisfy w_t − w_xx = 1 and w ≥ 0 on Q̄. The functions v± = ρw ∓ e then satisfy v±_t − v±_xx + c v± = ρ(1 + c w) ± r ≥ 0 and v± ≥ 0 on the parabolic boundary. The weak maximum principle for ∂_t − ∂_xx + c with bounded c ≥ 0 gives v± ≥ 0, so |e| ≤ ρw. Taking the smaller barrier pointwise finishes the proof. ∎

The checker verifies the data conditions exactly on the coefficient grid c[i][j] of x^i t^j:

  • row 0 vanishes (p(0, t) = 0);
  • every column sum vanishes (p(1, t) = 0);
  • column 0 equals the coefficients of u₀ (p(x, 0) = u₀).

The residual bound uses the tensor Bernstein enclosure of §7.1. Sensors must have 0 < x < 1 and 0 < t ≤ 1.

Corner compatibility. At a corner (c, 0) with c ∈ {0, 1}, p(c, t) = 0 forces p_t(c, 0) = 0, so the residual there is −u₀″(c) − g(c, 0). Any data with u₀″(c) + g(c, 0) ≠ 0 therefore force ρ to be at least that value. Polynomial candidates converge quickly only when the next order also holds. For g = 0, differentiating the equation in time at the lateral boundary gives u_xxt = 0. Consequently the equation differentiated twice in space requires u₀⁗(c) = (q f(u₀))″(c). When f′(0) = 0, this reduces to u₀⁗(c) = 0 (using u₀(c) = 0 and the odd-power form of f). A linear reaction term can make the right-hand side nonzero, so vanishing fourth derivatives alone need not suffice. The default cubic nonlinearity and u₀ = x − (5/3)x³ + x⁵ − (1/3)x⁶ satisfy both orders. In this setting, space-time Chebyshev candidates with n = 16 reach ρ ≈ 6.5·10⁻³, and n = 24 reaches ρ ≈ 4·10⁻⁴. The certified twins in tests/test_parabolic.py are θ₀ = 5 and θ₀ = 10 at noise level ε = 10⁻³.

7.3 Polynomial ODE systems

Let F(t, u; θ) be a polynomial vector field on ℝ × ℝⁿ and u₀(θ) a polynomial initial value. Signs are unrestricted. Consider

(O_θ)   u′ = F(t, u; θ),   u(0) = u₀(θ),   t ∈ [0, T],

with observables y = c·u + c₀ (affine in the state) measured at times tᵢ. F is smooth, so (O_θ) has a unique maximal solution, which may blow up. Write |·| for the max norm on ℝⁿ. Write μ(A) = maxᵢ (Aᵢᵢ + Σ_{j≠i} |Aᵢⱼ|) for its logarithmic norm, which is subadditive and positively homogeneous, hence convex.

A candidate consists of:

  • breaks 0 = b₀ < … < b_K ≥ T;
  • on each segment, polynomials p_k(τ) in τ = (t − b_k)/h_k;
  • a tube radius δ > 0.

On segment k, let ρ_k ≥ sup |p_k′ − F(t, p_k)|, and let μ_k ≥ sup { μ(∂F/∂u(t, v)) : t ∈ segment, |v − p_k(t)| ≤ δ }.

Theorem 11 (tube validation). Let E₀ = |p₀(0) − u₀|. For k ≥ 1, let E_k = Ē_{k−1} + |p_k(0) − p_{k−1}(1)|, where Ē_{k−1} is the end-of-segment bound below. Suppose that on every segment μ_k⁺ h_k < 1, where x⁺ = max(x, 0), and

(E_k + ρ_k h_k) / (1 − μ_k⁺ h_k) < δ.

Then the solution exists on [0, b_K], and for t = b_k + s in segment k

|u(t) − p_k(t)| ≤ E_k / (1 − μ_k s) + ρ_k s / (1 − μ_k⁺ s)  =: B_k(s),

with Ē_k = B_k(h_k).

Proof. Fix segment k and suppose |e(b_k)| ≤ E_k with e = u − p_k. As long as |e| ≤ δ, the segment from p_k to u lies in the tube. Hence e′ = A(t) e − r with A = ∫₀¹ ∂F/∂u(t, p_k + σe) dσ, and by convexity μ(A) ≤ μ_k. The upper Dini derivative then satisfies D⁺|e| ≤ μ(A)|e| + |r| ≤ μ_k|e| + ρ_k (Dahlquist 1958; see Söderlind, BIT 46 (2006) 631–652). Comparison with the scalar equation gives

|e(b_k + s)| ≤ E_k e^{μ_k s} + ρ_k (e^{μ_k s} − 1)/μ_k.

Two inequalities bound the right-hand side:

  • e^x ≤ 1/(1 − x) for every x < 1, including negative x.
  • (e^{μs} − 1)/μ = ∫₀ˢ e^{μσ} dσ ≤ s · e^{μ⁺s} ≤ s/(1 − μ⁺s).

Together they give B_k.

The supremum of B_k over [0, h_k] is at most (E_k + ρ_k h_k)/(1 − μ_k⁺ h_k) < δ. So the first time |e| reached δ would contradict the bound. The solution therefore stays in the tube, which also rules out blow-up, and the estimate holds on the whole segment. Induction over the segments, together with the triangle inequality at each jump, finishes the proof. ∎

Weighted norms. Given exact scales w_i > 0, apply the theorem to z_i = u_i / w_i. That system has right-hand side F_i / w_i and Jacobian J_ij w_j / w_i. In practice:

  • the residual and the initial or jump errors are divided by w_i;
  • the tube becomes |v_i − p_i| ≤ δ w_i;
  • μ_w = max_i (J_ii + Σ_{j≠i} |J_ij| w_j / w_i);
  • an observable c·u picks up Σ_j |c_j| w_j B.

Nothing else changes. Scaling matters for models whose states carry very different units (tests/test_ode.py::test_state_scales_rescue_badly_scaled_units).

Corollary. At a measurement (tᵢ, y = c·u + c₀, yᵢ) the checker compares |c·p(tᵢ) + c₀ − yᵢ| ± ‖c‖₁ B(tᵢ) with ε. The verdict rules are those of §7: a witness is VALID when every sum is ≤ ε, and INVALID when some difference exceeds ε. INVALID is sound because (O_θ) has a unique solution, which the theorem places within B of p.

Certified enclosures. Both ρ_k and the zero-order diagonal term of μ_k are univariate polynomials in τ, bounded by Bernstein enclosures with de Casteljau refinement (§6). The diagonal term uses a signed upper bound, so dissipative systems keep μ_k < 0.

Every other contribution is bounded through the tube expansion. Write v = p_k(τ) + δs with |s| ≤ 1 and expand ∂F/∂u(t, v) = Σ_α c_α(τ) s^α. Then |s^α| ≤ 1 gives

Jᵢᵢ + Σ_{j≠i} |Jᵢⱼ| ≤ sup c⁽ⁱⁱ⁾₀ + Σ_{α≠0} sup |c⁽ⁱⁱ⁾_α| + Σ_{j≠i} Σ_α sup |c⁽ⁱʲ⁾_α|.

Model strings are parsed by a closed grammar (twindex/mpoly.py): numbers, names, + − *, division by nonzero constants, and integer powers. Decimal literals are read exactly, so "0.1" means 1/10.

Example (tests/test_ode.py). Consider SIR (S′ = −βSI, I′ = βSI − γI, S₀ = 0.99, I₀ = 0.01), with I observed at t = 0.5, 1, …, 4 and ε = 10⁻³. Twindex certifies both (β, γ) = (0.5, 0.2) and (0.72, 0.4178) as data-compatible. These share the early growth rate βS₀ − γ, yet their reproduction numbers β/γ are 2.5 and about 1.72. Early-epidemic data of this kind provably cannot determine R₀.

7.4 Elementary right-hand sides (certified polynomialization)

Right-hand sides may use + − *, division by any expression, integer powers (including negative ones), exp, sin and cos. This covers Michaelis–Menten and Hill kinetics, Monod growth, Holling responses, frequency-dependent incidence βSI/N, Arrhenius rates exp(−E/T) and periodic forcing.

Construction. Atoms are processed innermost first. Each non-polynomial atom a whose argument depends on time or on a state becomes an auxiliary state. Write ḋ = ∂_t d + Σ_x ∂_x d · F̃_x for the derivative of an argument d along the lifted vector field F̃; the sum runs over every lifted state x. Then:

atom auxiliary state equation initial value
1/D z z′ = −z² Ḋ 1/D(0) (exact rational)
exp(g) z z′ = z ġ exp(g(0)) (exact rational interval)
sin(g), cos(g) s, c s′ = c ġ, c′ = −s ġ sin g(0), cos g(0) (intervals)

A reciprocal 1/D whose argument involves only parameters is a derived constant 1/D(θ). It is computed exactly per θ, and its value is INVALID if D(θ) = 0. The lifted right-hand side F̃ is a polynomial.

Theorem 12. Let v = (u, z) solve the lifted system on [0, T] with the exact (possibly irrational) initial values above. Then u solves the original system on [0, T], and every denominator D(u(t)) is nonzero there.

Proof. Take the atoms in creation order. Each atom's argument involves only states and earlier atoms, which are already known to equal their intended expressions.

  • Reciprocal. Let w = D z. Then w′ = Ḋ z + D(−z² Ḋ) = Ḋ z (1 − w), and w(0) = 1. This is a linear equation in w with continuous coefficient, so w ≡ 1. Hence z = 1/D and D ≠ 0.
  • Exponential. w = z e^{−g} has w′ = z ġ e^{−g} − z e^{−g} ġ = 0, so w ≡ 1.
  • Sine and cosine. Let φ = g(u(t)), so that φ′ = ġ. Both (s, c) and (sin φ, cos φ) solve the linear system s′ = c ġ, c′ = −s ġ with the same initial value. By uniqueness they coincide, and the phase is pinned, not only the circle s² + c² = 1.

So each atom equals its intended expression, F̃(v) = F(u), and u′ = F(u). ∎

Certification. Theorem 11 applies to the lifted polynomial system with one change: an initial value given only as an interval [lo, hi] contributes max(|p(0) − lo|, |p(0) − hi|) to E₀. This bounds |p(0) − v(0)| for the true initial value whatever it is. Nothing is chosen inside the interval, so a candidate that assumes a wrong value of, say, exp(1) pays for the difference.

The enclosures of exp, sin and cos at a rational point x with |x| ≤ 64 are exact Taylor sums with Lagrange remainder |x|^k/k!, multiplied by 3^{⌈x⌉} for exp with x > 0. They are rounded outward to multiples of 2⁻¹⁶⁰. Initial arguments that involve earlier atoms are evaluated in exact interval arithmetic.

The checker derives the lifting itself, deterministically, from the original model strings in the certificate. The witness supplies candidates for the auxiliary states as well. Observables refer to the user's states only. A model whose true solution reaches a singularity, for example u′ = 1/(1 − u) at t = 1/2, can never be certified past it, because Theorem 11 would prove that the lifted solution exists there.

Prior work. Polynomialization of elementary ODEs is classical; see Hemery, Fages and Soliman (2021) and the quadratization literature. What is new here is that the lifting runs inside an exact, replayable certificate checker.

8. What a certificate does not say

  • It concerns the declared data y and the bounded-noise model max_i |noise_i| ≤ ε. It is not a statistical statement and does not assert non-identifiability for exact (noise-free) data.
  • A pair certificate gives a lower bound on D_Q(y). Upper bounds need a covering of Θ by a paving certificate (Theorem 5), and a paving is only as tight as its unresolved leaves allow.
  • INCONCLUSIVE means the sufficient inequalities did not close. It is not evidence that the claim is false.

The ingredients are classical; Twindex's contribution is putting them together in an exact, independently replayable checker.

  • Optimal recovery and the radius of information: Micchelli & Rivlin (1977); Foucart & Liao, arXiv 2111.02601 and arXiv 2401.11112.
  • Optimal recovery for parametric PDEs: Binev et al., SIAM/ASA JUQ 5 (2017); Cohen, Dahmen, Mula & Nichols, SIAM/ASA JUQ 10 (2022), arXiv 2009.02687.
  • Guaranteed (set-membership) parameter estimation and set inversion: Milanese & Vicino (1991); Jaulin & Walter, Automatica 29 (1993); Raïssi, Ramdani & Candau, Automatica (2004).
  • Certified reduced-basis possibility regions: Sen et al., J. Comput. Phys. (2006); Nguyen, Rozza, Huynh & Patera (2011).
  • Maximum-principle a posteriori analysis for semilinear problems: Nochetto, Schmidt, Siebert & Veeser.
  • Computer-assisted proofs for PDEs: Nakao, Plum & Watanabe (2019).
  • Error certification for physics-informed networks: Eiras et al., ICML 2024.