Certified Baker–Davenport reduction for the unit right-hand-side case

[roed-math]

Goal

For matched simultaneous-Pell systems in the epic's supported class (parent #58, "Shared mathematical contract") with unit right-hand sides — (N1, N2) in {1, -1}^2 and both seed sets nonempty — replace the bounded orbit search by a proof: parametrize each equation's solutions by finitely many one-sided rays (constructed and proved complete below), derive an explicit upper bound on the ray indices from a linear form in three logarithms (Matveev's theorem in the form of [BugeaudMignotteSiksek2006, Thm 9.4], transcribed in full below), reduce it with the Dujella–Pethő lemma ([DujellaPetho1998, Lemma 5a], stated in full below and implemented in this child on #72's certified-arithmetic utilities), finish with the existing orbit walk below the final bound, and return the full finite list of all signed integer solutions with completeness="proved".

Target base

• Repository state: branch review-architecture @ 85f8e0a (the composed
  tree on the roed-math fork: wave-1 code = PRs Packaging, CI, design notes, and the annotated bibliography #2–Family: waring — Waring-type diagonal representation #48, plus the wave-2 and
  architecture-review series, which are not yet opened as PRs).
• Required merged PRs: Packaging, CI, design notes, and the annotated bibliography #2–Family: waring — Waring-type diagonal representation #48 (wave 1); the wave-2/architecture-review series
  must land after them before this issue's work can merge to main.
  Implementation happens on the composed tree, not on a split branch.

Dependencies

• Parent: Rigorous completion of simultaneous Pell searches #58 (shared mathematical contract binding).
• Blocked by: Complete seed-orbit enumeration for generalized Pell equations #80 (seed-orbit enumeration) and Rigorous de Weger completion of the exponential searches #55's Certified logarithmic-bound and LLL-reduction primitives #72 (certified ball
  evaluation, precision escalation, and the shared round-record/certificate
  conventions — shared infrastructure only, separate mathematical
  reduction: the S-unit/de Weger pipeline and this Baker–Davenport-style
  reduction are different mathematics. This child does not call Certified logarithmic-bound and LLL-reduction primitives #72's
  generic reduce_linear_form, and nothing here may describe Certified logarithmic-bound and LLL-reduction primitives #72's de Weger
  lemmas as the Dujella–Pethő lemma: the Dujella–Pethő inequality is
  implemented in this child, using Certified logarithmic-bound and LLL-reduction primitives #72 only for rigorous ball evaluation of
  exact values, the precision-escalation discipline, and the unified round-
  record field names shared with Evidence and certificate integration for simultaneous Pell bounds #82).
• Unblocks: Evidence and certificate integration for simultaneous Pell bounds #82, Extension to supported nonunit right-hand sides #83, Simultaneous-Pell solver integration and partial fallback boundaries #84.

Supported inputs

The epic contract restricted to the unit regime:

• shared unknown = the B-coordinate of both equations (x^2 - D1*z^2 = N1; y^2 - D2*z^2 = N2 up to names/orientation);
• D1, D2 concrete positive nonsquares (squarefree not required), with QQ(sqrt(D1)) != QQ(sqrt(D2)) — exact integer test: D1*D2 is not a perfect square; the dependent case falls back (Simultaneous-Pell solver integration and partial fallback boundaries #84);
• (N1, N2) in {1, -1}^2;
• both seed sets nonempty (Complete seed-orbit enumeration for generalized Pell equations #80): for N = 1 the seed set always contains (1, 0); for N = -1 it may be empty — if the first equation's seed set is empty the existing proved-empty branch answers (unchanged, Simultaneous-Pell solver integration and partial fallback boundaries #84), and if only the second one's is, this child makes no claim and the system stays on the fallback path (Simultaneous-Pell solver integration and partial fallback boundaries #84 lists it in the boundary).

Everything outside this regime is untouched by this child: the dispatch must fall through to the current code path byte-for-byte.

Output and API contract

• solve() on a system in the regime returns a SolutionSet with:
  • kind="finite-complete", completeness="proved" (so the derived complete view is True), scope="all integer solutions";
  • solutions = the full list of signed integer solution tuples, ordered by the unknowns' first appearance in the input and sorted ascending (the existing sorted(...) convention);
  • every point passing the solvers._verified firewall against the original ParsedSystem (the firewall stays in place, never bypassed);
  • typed evidence (this child, provisional in structure, final in mathematical content): a kind="theorem" record for the initial bound naming Matveev via BMS Theorem 9.4 with reference="BugeaudMignotteSiksek2006", a kind="theorem" record for the reduction lemma with reference="DujellaPetho1998", and a kind="external-computation" record carrying the run's numeric bound chain. Evidence and certificate integration for simultaneous Pell bounds #82 replaces the payload structure with the shared certificate schema; do not invent a competing schema here beyond what Evidence and certificate integration for simultaneous Pell bounds #82 specifies. Reduction-round entries use exactly the unified Certified logarithmic-bound and LLL-reduction primitives #72/Evidence and certificate integration for simultaneous Pell bounds #82 field names {"dimension", "scaling_C", "precision_bits", "b1_norm_lower", "new_bound"} with the Dujella–Pethő semantics pinned under "Reduction and termination" below.
• New per-equation primitive (used by this child and Extension to supported nonunit right-hand sides #83): one_sided_rays(orbits) added to diophantine_classifier/pell_orbits.py, taking a Complete seed-orbit enumeration for generalized Pell equations #80 SeedOrbits record and returning a tuple of frozen Ray records (seed, sign, shift, base) — exact ZZ data, one ray per seed, in seed order — implementing the deterministic construction below. It changes nothing in Complete seed-orbit enumeration for generalized Pell equations #80's existing contract.
• The docstring EXAMPLES of _solve_simultaneous_pell are updated: the existing doctest asserting S.complete is False for x^2 - 3*z^2 = 1; y^2 - 8*z^2 = 1 flips to the new proved behavior (True, with the full list).
• CAPABILITIES["simultaneous-pell"] must declare completeness {"partial", "proved"} and kinds {"witness", "finite-complete", "empty"}. At the target base it already declares exactly this — the "proved" value was until now exercised only by the first-equation-empty branch; after this child it is also exercised by finite-complete pipeline results. The entry itself (domains, goals, kinds, completeness) stays unchanged.

Algorithm and conventions

Orbit coordinates (#80)

For equation i write eps_i = u_i + v_i*sqrt(D_i) > 1 for the fundamental automorph (fundamental solution of u^2 - D_i*v^2 = 1; N(eps_i) = 1, so eps_i^(-1) = eps_i', the algebraic conjugate) and identify a solution (x, z) of x^2 - D*z^2 = N with the element beta = x + z*sqrt(D) of ZZ[sqrt(D)] (beta' = x - z*sqrt(D), beta*beta' = N; the identification is a bijection because sqrt(D) is irrational). #80 supplies the finite normalized seed list; the full solution set is the disjoint union over seeds s of the classes { ±eps^k * (x_s + z_s*sqrt(D)) : k in ZZ }. Everything below is per equation until "Smallness at shared values".

One-sided orbit parametrization

This subsection converts the two-sided k in ZZ classes into finitely many one-sided rays beta_r * eps^n (n >= 0) and proves the exactly-one coverage that #81/#83's per-ray analysis requires. All numeric claims in this subsection are machine-checked (verified against Sage 10.7 / PARI 2.17.2; fixtures below).

From two-sided to one-sided (the conjugation calculation). Fix a seed s = (x_s, z_s) with alpha = x_s + z_s*sqrt(D). Because N(eps) = eps * eps' = 1, we have eps^(-1) = eps', hence for every k >= 0

alpha * eps^(-k) = alpha * (eps')^k = ((alpha)' * eps^k)'   —   i.e.  alpha*eps^(-k) = (alpha'*eps^k)'.

Therefore the two-sided class splits into two one-sided sequences plus signs:

{ ±alpha*eps^k : k in ZZ }  =  { ±alpha*eps^n : n >= 0 }  ∪  { ±(alpha'*eps^n)' : n >= 1 }.

Conjugation acts on coordinates as (x, z)' = (x, -z) and negation as -(x, z) = (-x, -z); both lie in the coordinate sign group (x, z) -> (±x, ±z) under which the reported solution set is closed (parent contract: full signed tuples). So, up to the allowed sign symmetries, the class's tuples are exactly the tuples of the two one-sided sequences alpha*eps^n (n >= 0) and alpha'*eps^n (n >= 1). The subset of {±alpha, ±alpha'} * eps^(n>=0) that this uses is: one sign choice of alpha for the forward tail, and one sign choice of alpha' for the backward tail (shifted to n >= 1); the remaining sign choices are recovered by the sign group. Listing these raw sequences directly would, however, double-count small elements between a ray and its conjugate ray (concrete repeat, D = 13, N = 27: the seed (-12, 3)'s sequence at n = 0 sign-normalizes to (12, 3), which is the n = 0 element of the different seed (12, 3)'s ray — verified against Sage 10.7). The anchored construction below eliminates every such repeat and is the parametrization this epic pins.

Quadrant normalization. Let Q := {(x, z) : x >= 0, z > 0}. Every solution with z != 0 has exactly one image in Q under the sign group (the four images are (±|x|, ±|z|); exactly (|x|, |z|) lies in Q), and z = 0 solutions have none (handled separately below). For a shared value analysis this is the normalization z > 0, x_i = |x_i| of the parent contract; note beta = x + z*sqrt(D) >= sqrt(D) > 1 on Q.

Deterministic construction (per seed, exact integer arithmetic only). Let E(x, z) = (u0*x + D*v0*z, v0*x + u0*z) (multiplication by eps). For each #80 seed s = (x_s, z_s), in the order of SeedOrbits.seeds:

1. sigma := the sign of x_s + z_s*sqrt(D), decided exactly: if x_s >= 0 and z_s >= 0 (not both 0): +1; if x_s <= 0, z_s <= 0: -1; if x_s > 0, z_s < 0: sign(N); if x_s < 0, z_s > 0: -sign(N). (Proof of the mixed cases: for x > 0 > z, x + z*sqrt(D) > 0 iff x^2 > D*z^2 iff N > 0; symmetrically for the other.)
2. (x, z) := sigma * (x_s, z_s); t := 0. While not (x >= 0 and z > 0): (x, z) := E(x, z); t := t + 1.
3. Ray(seed=s, sign=sigma, shift=t, base=(x, z)), i.e. beta_r := sigma * (x_s + z_s*sqrt(D)) * eps^t, the ray being { beta_r * eps^n : n >= 0 }.

The loop terminates and the result is exactly the minimal-z element of Q ∩ class(s) (Lemma R2); moreover t = 0 if and only if the seed itself lies in Q (proof below), so t is 0 or the small positive count of forward steps — never negative, and no backward walk is needed.

Lemma R1 (positive part of a class is one <eps>-orbit). eps^k > 0 for all k, so the two orbits {alpha*eps^k} and {-alpha*eps^k} consist of elements of constant sign, opposite to each other. The positive elements of the class are { a*eps^k : k in ZZ } with a := sigma*alpha > 0 as in step 1–2. Only positive elements can lie in Q (beta > 0 there).

Lemma R2 (the quadrant cuts a single one-sided tail). Write x_k + z_k*sqrt(D) = a*eps^k (k in ZZ), so x_k = (a*eps^k + a'*eps^(-k))/2 and z_k = (a*eps^k - a'*eps^(-k))/(2*sqrt(D)), with a*a' = N.

• If N > 0: a' > 0, hence x_k > 0 for all k, and (viewing k as a real variable) d/dk[ a*eps^k - a'*eps^(-k) ] = log(eps) * (a*eps^k + a'*eps^(-k)) = 2*log(eps)*x_k > 0, so z_k is strictly increasing in k, tending to -infinity and +infinity. Therefore {k : (x_k, z_k) in Q} = {k : z_k > 0} = {k >= n_min} for a unique n_min in ZZ; algebraically z_k > 0 iff eps^(2k) > a'/a. (z_k = 0 can occur for at most one k, exactly when a'/a is an even power of eps; then N = x_k^2 is a perfect square and that element is a z = 0 solution, outside Q.)
• If N < 0: a' < 0, hence z_k > 0 for all k, and d(x_k)/dk = log(eps)*(a*eps^k - a'*eps^(-k))/2 = log(eps)*sqrt(D)*z_k > 0, so x_k is strictly increasing from -infinity to +infinity. Therefore Q ∩ class = {k : x_k >= 0} = {k >= n_min}; algebraically x_k >= 0 iff eps^(2k) >= -a'/a. (x_k = 0 occurs for at most one k, exactly when D | N with |N|/D a perfect square; such an element lies in Q and is a legitimate ray element.)

In both cases Q ∩ class(s) is the single one-sided tail { a*eps^k : k >= n_min }; re-indexing by n = k - n_min gives the ray { beta_r * eps^n : n >= 0 } with beta_r = a*eps^(n_min). The construction's walk finds it: if sigma*s in Q then, because the seed has the minimal |z| of its whole class (#80's rule) and z determines the Q-element uniquely (x = +sqrt(N + D*z^2)), sigma*s is already the minimal-z element of Q ∩ class — the base, with t = 0; otherwise sigma*s sits strictly below the tail and the forward walk crosses into Q at exactly the tail's first element (monotone single crossing, by the monotonicity just proved), so the base is again the minimal-z Q-element.

Lemma R3 (exactly-one coverage). Every positive z-value taken by any solution of x^2 - D*z^2 = N occurs in exactly one ray at exactly one index. Proof: given such a z, the unique Q-solution with that z is q = (+sqrt(N + D*z^2), z). The classes partition the solutions, so q lies in exactly one class; by Lemma R2 that class's Q-part is its single ray, and membership at two indices is impossible since n -> beta_r*eps^n is strictly increasing (eps > 1). Distinct rays therefore have disjoint z-sets, and no deduplication rule is needed. The full signed solution set is recovered as { (±x, ±z) : (x, z) a ray element } ∪ { (±sqrt(N), 0) if N is a perfect square }. (What became of the backward tail: the elements a*eps^k with k < n_min are not in Q; their sign-normalized images lie — by the conjugation calculation above — in the ray of the conjugate class, the class of ±alpha', which is another seed's class (or the same, for a self-conjugate class); the anchored rays list each such image exactly once there. This is precisely why the naive {alpha, alpha'}*eps^(n>=0) lists would repeat small-n elements while the anchored rays do not.)

Binet form, monotonicity index, positivity. On ray r, for n >= 0:

z_n = (beta_r * eps^n - beta_r' * eps^(-n)) / (2*sqrt(D)),      x_n = (beta_r * eps^n + beta_r' * eps^(-n)) / 2,

with beta_r * beta_r' = N, beta_r >= sqrt(D) > 1, and Lambda := beta_r * eps^n > 0 — the positivity hypothesis every analytic step below uses. z_n is strictly increasing from n0(r) = 0 on: z_{n+1} = v0*x_n + u0*z_n >= u0*z_n >= 2*z_n > z_n since x_n >= 0, z_n > 0, u0 >= 2 on Q — so the explicit monotonicity index is n0(r) = 0 for every ray, by the anchoring (this is a proved property of the construction, not an assumption). Also x_{n+1} = u0*x_n + D*v0*z_n > x_n, and beta_r*eps^n >= max(x_n, z_n) bounds both coordinates.

The N > 0, N < 0, and z = 0 cases, summarized. N > 0: all ray elements have x > 0; z = 0 solutions exist iff N is a perfect square, are exactly (±sqrt(N), 0), lie in no ray, and are reported separately (in the system: found by the finishing walk, which streams them with the #80 seeds). N < 0: all class elements on the positive side already have z > 0; a ray may contain one x = 0 element (iff D | N, |N|/D square); z = 0 is impossible. In this child's unit regime (N in {1, -1}): for N = 1 the single seed is (1, 0) and the ray is {eps^n : n >= 1} with z = 0 giving the extra tuples (±1, 0); for N = -1 there is at most one class and no z = 0 or x = 0 element (x^2 - D*z^2 = -1 forces x*z != 0).

Deterministic ray order. Rays are listed in the order of their seeds in SeedOrbits.seeds (#80's pinned (|y|, y, x) sort). Ray pairs for the system are enumerated in lexicographic seed-pair order, exactly as the certificate (#82/#83) records them.

Verified fixtures (both machine-checked end-to-end; verified against Sage 10.7 / PARI 2.17.2). For D = 13, N = 27 (fundamental (649, 180), four #80 seeds):

seed (-12, 3)  -> sign -1, shift 1, base (768, 213)     z-values 213, 276477, ...
seed (12, 3)   -> sign +1, shift 0, base (12, 3)        z-values 3, 4107, ...
seed (-40, 11) -> sign -1, shift 1, base (220, 61)      z-values 61, 79189, ...
seed (40, 11)  -> sign +1, shift 0, base (40, 11)       z-values 11, 14339, ...

Enumerating all signed solutions with |z| <= 10^6 by brute force (8 sign-normalized representatives, z in {3, 11, 61, 213, 4107, 14339, 79189, 276477}) and matching each against the rays confirms: every representative occurs in exactly one ray at exactly one index, the four ray z-sets are pairwise disjoint, every listed ray element solves the equation with x >= 0 < z, and there is no z = 0 solution (27 is not a square). For the unit case D = 3, N = 1 (fundamental (2, 1), single seed (1, 0)): sign +1, shift 1, base (2, 1); the ray's z-values 1, 4, 15, 56, 209, 780, 2911, 10864, 40545, 151316, 564719 cover all eleven sign-normalized representatives with 0 < |z| <= 10^6 exactly once, and the z = 0 solutions are exactly (±1, 0). The identity alpha*eps^(-k) = (alpha'*eps^k)' was checked exactly in QuadraticField(D) for both fixtures.

Smallness at shared values (per ray pair)

Fix rays r1 (equation 1) and r2 (equation 2) with bases beta_1, beta_2. A shared value z > 0 means z = z_j on r1 and z = z_k on r2 for unique j, k >= 0 (Lemma R3); write Lambda_1 = beta_1*eps_1^j, Lambda_2 = beta_2*eps_2^k (both > 1). From Lambda_i - N_i/Lambda_i = 2*z*sqrt(D_i) (both equal because the z's match):

sqrt(D2)*Lambda_1 - sqrt(D1)*Lambda_2 = sqrt(D2)*N1/Lambda_1 - sqrt(D1)*N2/Lambda_2,

so with t := sqrt(D2)*Lambda_1/(sqrt(D1)*Lambda_2) - 1 = (sqrt(D2)*N1/Lambda_1 - sqrt(D1)*N2/Lambda_2) / (sqrt(D1)*Lambda_2) and the quadrant bounds Lambda_1 >= z*sqrt(D1), Lambda_2 >= z*sqrt(D2) (valid since x_i >= 0):

|t| <= sqrt(D2)*|N1|/(sqrt(D1)*Lambda_1*Lambda_2) + |N2|/Lambda_2^2 <= |N1|/(D1*z^2) + |N2|/(D2*z^2) = T0/z^2,
T0 := |N1|/D1 + |N2|/D2      (exact rational),

with no side condition — the inequality holds at every shared z >= 1. In this unit regime T0 = 1/D1 + 1/D2. Checked numerically at the shared values of all committed fixtures, e.g. (3,1)x(8,1), z=1: |t| = 0.0456 <= 11/24; (13,27)x(2,7), z=3: |t| = 0.0317 <= 0.6197 (verified against Sage 10.7).

Two exact consequences (both verified on the fixture rays, Sage 10.7):

• Index decay. Let j0(r1) := the least j >= 0 with beta_1^2 * eps_1^j >= 2*|N1| (an exact comparison in ZZ[sqrt(D1)]; j0 = 0 on every committed fixture ray). For j >= j0, |beta_1'*eps_1^(-j)| = |N1|/(beta_1*eps_1^j) <= beta_1*eps_1^j/2, so z_j >= beta_1*eps_1^j/(4*sqrt(D1)) and hence

  |t| <= (16*T0*D1/beta_1^2) * eps_1^(-2j)          for j >= j0(r1).
  

• The three-logarithm forms. With gamma_1 := eps_1, gamma_2 := eps_2, gamma_3 := beta_1*sqrt(D2)/(beta_2*sqrt(D1)) (the ray constant) and (b_1, b_2, b_3) = (j, -k, 1):

  gamma_1^j * gamma_2^(-k) * gamma_3 - 1 = t          (exactly — this is Matveev's Lambda),
  log(gamma_1^j*gamma_2^(-k)*gamma_3) = j*log(eps_1) - k*log(eps_2) + log(gamma_3) = log(1 + t),
  

  so the alpha^b - 1 form the theorem below bounds is t; no exp/log conversion enters the Matveev step. For the Dujella–Pethő step divide the log form by log(eps_2): whenever z^2 >= 2*T0 (so |t| <= 1/2 and |log(1+t)| <= 2|t|),

  |j*theta - k + mu| <= (2/log(eps_2)) * |t| <= A * B^(-j),
  theta := log(eps_1)/log(eps_2),   mu := log(gamma_3)/log(eps_2),
  A := 32*T0*D1 / (beta_1^2 * log(eps_2)),   B := eps_1^2,       valid for j >= max(j0(r1), jz(r1)),
  

  where jz(r1) := the least index with z_{jz}^2 > 2*T0 on r1 (exact integer test on the actual ray elements). The finitely many indices below max(j0, jz) are absorbed into the reduction bounds below and covered by the finishing walk.

Nonvanishing, decided exactly (no implementation-time deferral). t = 0 iff sqrt(D2)*Lambda_1 = sqrt(D1)*Lambda_2. Writing Lambda_1 = P + Q*sqrt(D1), Lambda_2 = R + S*sqrt(D2) (P, Q, R, S in ZZ) and expanding in the basis {1, sqrt(D1), sqrt(D2), sqrt(D1*D2)} of the degree-4 field QQ(sqrt(D1), sqrt(D2)) (degree 4 exactly because D1*D2 is not a square — the regime guarantee; degree checked in Sage for the fixture pairs), the identity forces P = R = 0 and Q = S: that is, t = 0 happens exactly at a shared value with x_1 = x_2 = 0. Such a point requires N_i = -D_i * (square) for both equations — impossible in this unit regime (x_i = 0 forces N_i = -D_i*z_i^2 <= -2, excluded by |N_i| = 1) — and in the general regime (#83) it is a single explicit index pair with z^2 = |N1|/D1 <= T0 < 2*T0, hence below the jz floor and inside the finishing walk. The lower-bound theorem is therefore only ever applied to a form proven nonzero, and the Dujella–Pethő inequality's strict 0 < ... hypothesis holds on the indices it is applied to.

Initial bound — Matveev's theorem in the form of [BugeaudMignotteSiksek2006, Thm 9.4]

Statement transcribed with the source open (ar5iv rendering of arXiv:math/0403046, Section 9.1; notation as in the source):

Let L be a number field of degree D, let alpha_1, ..., alpha_n be non-zero elements of L and b_1, ..., b_n be rational integers. Set

B = max{ |b_1|, ..., |b_n| },        Lambda = alpha_1^(b_1) * ... * alpha_n^(b_n) - 1.

Let h denote the absolute logarithmic height and let A_1, ..., A_n be real numbers with

A_j >= h'(alpha_j) := max{ D*h(alpha_j), |log alpha_j|, 0.16 },     1 <= j <= n.

We call h' the modified height (with respect to the field L).

Theorem 9.4. Assume that Lambda is non-zero. We then have

log |Lambda| > -3 * 30^(n+4) * (n+1)^(5.5) * D^2 * (1 + log D) * (1 + log(n*B)) * A_1 * ... * A_n.

Furthermore, if L is real, we have

log |Lambda| > -1.4 * 30^(n+3) * n^(4.5) * D^2 * (1 + log D) * (1 + log B) * A_1 * ... * A_n.

Specialization arithmetic (n = 3, real field). Per ray pair take L = QQ(sqrt(D1), sqrt(D2)) (a real field), alpha_1 = eps_1, alpha_2 = eps_2, alpha_3 = gamma_3, (b_1, b_2, b_3) = (j, -k, 1), so Lambda = t != 0 (nonvanishing above) and B = max(j, k, 1). The degree is D = 4: D1*D2 not a square is exactly the nondegeneracy condition making [QQ(sqrt(D1), sqrt(D2)) : QQ] = 4 (checked in Sage for (3,8), (13,2), (2,3)). The real-field bound applies with the explicit prefactor

1.4 * 30^6 * 3^(4.5) * 4^2 * (1 + log 4) = 5.46695133402655...e12    (certified ball value; upper bound 5.4669513341e12, verified against Sage 10.7),

so log |t| > -kappa_M * (1 + log B) with kappa_M := 5.46695...e12 * A_1 * A_2 * A_3, all factors evaluated as certified upper bounds at run time (an upper bound for each A_j keeps the theorem valid, since the hypothesis is one-sided: A_j >= h'(alpha_j)).

Exact heights for the A_j (all formulas derived below and machine-checked against Sage's exact global_height on the fixture data — verified against Sage 10.7).

• Quadratic units: eps_i has minimal polynomial T^2 - 2*u_i*T + 1, conjugate 1/eps_i in (0, 1), so its Mahler measure is eps_i and h(eps_i) = log(eps_i)/2 exactly. Hence D*h(eps_i) = 4 * log(eps_i)/2 = 2*log(eps_i) >= |log eps_i| and >= 0.16 (every fundamental automorph has u_i >= 2, v_i >= 1, D_i >= 2, so eps_i >= 2 + sqrt(2) and 2*log(eps_i) > 2.45), so take A_1 = 2*log(eps_1), A_2 = 2*log(eps_2) (certified upper balls). Checked: h((2+sqrt3)) = 0.658478948... = log(2+sqrt3)/2, likewise for (3,1|D=8), (649,180|D=13), (3,2|D=2).

• Ray bases: beta = x_r + z_r*sqrt(D) is an algebraic integer with minimal polynomial T^2 - 2*x_r*T + N, so h(beta) = (log^+|beta| + log^+|beta'|)/2 = (log(beta) + log^+(|N|/beta))/2 exactly (beta > 1 on Q). Checked on five fixture bases (e.g. h(12 + 3*sqrt13) = 1.6479184330).

• h(sqrt(D)) = log(D)/2 exactly (minimal polynomial T^2 - D, Mahler measure D).

• The ray constant gamma_3 = beta_1*sqrt(D2)/(beta_2*sqrt(D1)): exact value by the Mahler formula — compute the minimal polynomial of gamma_3 over QQ exactly in the biquadratic field, clear denominators to a primitive ZZ polynomial P, then

  h(gamma_3) = ( log|lc(P)| + sum over roots rho of P of log^+|rho| ) / deg(P)
  

  with certified root enclosures (this is the standard Mahler-measure height formula; checked to agree with Sage's global_height for both fixture gamma_3's, e.g. for the classical fixture gamma_3 = (2+sqrt3)*sqrt8/((3+sqrt8)*sqrt3): primitive polynomial 9*T^4 + 288*T^3 - 1104*T^2 + 768*T + 64, h = 1.687043381954, |log gamma_3| = 0.0446253...). A cheaper valid alternative the implementation may use for A_3 is the subadditive upper bound h(gamma_3) <= h(beta_1) + h(beta_2) + log(D1)/2 + log(D2)/2 (from h(xy) <= h(x)+h(y), h(1/x) = h(x); for the classical fixture it evaluates to 3.1288794507 >= 1.6870... — checked). Either way A_3 := max(4*h_upper(gamma_3), |log gamma_3|_upper, 0.16).

Closing the inequality. Combining log|t| > -kappa_M*(1 + log B) with the index decay log|t| <= log(16*T0*D1/beta_1^2) - 2*j*log(eps_1) (valid for j >= j0), and eliminating k via k <= theta*j + (|log gamma_3| + log 2)/log(eps_2) (valid for j >= jz; from the log form with |log(1+t)| <= log 2), so B <= max(1, j, theta_up*j + c_up):

2*j*log(eps_1) <= log(16*T0*D1/beta_1^2) + kappa_M * (1 + log(max(1, j, theta_up*j + c_up))),

an inequality of the shape x <= a*log(x) + b resolved by certified doubling-plus-bisection (#72 ball comparisons, upper roundings) into

B0(r1, r2) := max( the least bound beyond which the inequality certifiably fails,  j0(r1),  jz(r1) )

— an explicit integer with: every shared value of the pair has ray index j <= B0. (Illustrative magnitude for the classical fixture: kappa_M = 3.4257e14 and B0 ≈ 4.83e15 — verified against Sage 10.7; the committed exact value is locked by the #82 certificate from the first verified run.) The k-range for the certificate box is k <= floor(theta_up*B0 + c_up).

Reduction and termination — the Dujella–Pethő lemma, implemented here

Statement, transcribed with the source open (A. Dujella and A. Pethő, A generalization of a theorem of Baker and Davenport, Quart. J. Math. Oxford (2) 49 (1998), 291–306, Section 5, Lemma 5 — transcription from the author-hosted preprint at https://web.math.pmf.unizg.hr/~duje/pdf/qjm.pdf, cross-checked verbatim against the restatement in arXiv:1510.05579, Lemma 11, which cites it as "[6, Lemma 5a]"):

Lemma 5. Suppose that M is a positive integer. Let p/q be the convergent of the continued fraction expansion of kappa such that q > 6M and let eps = ||mu*q|| - M*||kappa*q||, where || . || denotes the distance from the nearest integer.

a) If eps > 0, then there is no solution of the inequality

0 < m*kappa - n + mu < A*B^(-m)                                   (16)

in integers m and n with

log(A*q/eps)/log(B) <= m <= M.

b) Let r = floor(mu*q + 1/2). If p - q + r = 0, then there is no solution of inequality (16) in integers m and n with max(log(3*A*q)/log(B), 1) < m <= M.

Only part a) is used by this pipeline. Proof sketch (so the citation is attribution, not a load-bearing pointer; this follows the original's half-page proof): multiply 0 < m*kappa - n + mu < A*B^(-m) by q and write m*kappa*q = m*(kappa*q - p) + m*p; then mu*q - (n*q - m*p) differs from the (positive, < q*A*B^(-m)) left side by m*(kappa*q - p), so ||mu*q|| <= |mu*q - (n*q - m*p)| < q*A*B^(-m) + m*|kappa*q - p|, because n*q - m*p is an integer. For a convergent p/q, |kappa*q - p| = ||kappa*q||; with m <= M this gives q*A*B^(-m) > ||mu*q|| - M*||kappa*q|| = eps, and if eps > 0, taking logarithms yields m < log(A*q/eps)/log(B). ∎ Two-sided corollary (needed because our form has no fixed sign): under the same hypotheses there is no solution of 0 < |m*kappa - n + mu| < A*B^(-m) with log(A*q/eps)/log(B) <= m <= M — for the negative-sign case run the same three-step estimate with (-kappa, -mu, -p, -n) in place of (kappa, mu, p, n); every quantity ||mu*q||, ||kappa*q||, and hence eps, is unchanged under negation. ∎

Application data (per ray pair; all values from "Smallness at shared values"). kappa = theta = log(eps_1)/log(eps_2) — irrational in this regime (if theta = p/q then eps_1^q = eps_2^p would be a rational integer unit > 1 in QQ(sqrt(D1)) ∩ QQ(sqrt(D2)) = QQ, impossible); mu = log(gamma_3)/log(eps_2); A = 32*T0*D1/(beta_1^2*log(eps_2)) (certified upper); B = eps_1^2 (certified lower for log B); m = j, n = k, M = the current bound B_r >= j. The strict positivity 0 < |j*theta - k + mu| holds on every index the corollary is applied to (nonvanishing above). Conclusion of one round:

B_{r+1} := max( floor( upper( log(A*q/eps_lo) / log(B) ) ),  j0(r1),  jz(r1) )        with certified eps_lo <= eps, eps_lo > 0.

Certified mechanics (this is where #72 is used — and all that it is used for): theta, mu, A, log B are #72 exact-value/ball evaluations at the escalating working precision; the continued-fraction convergents of theta are read off the ball enclosure (partial quotients emitted only while both endpoints of the enclosure agree on them — otherwise precision doubles); the round takes the first convergent with q > 6*M whose eps_lo (an exact rational lower bound extracted from the balls of ||mu*q|| and ||theta*q||) is certified > 0, else the next convergent, under #72's max_precision_bits cap — the exact convergent-hunting practice of the source paper ("if the first convergent such that q > 6M does not satisfy the condition ... then we use the next convergent"). Cap reached without a certified round = pipeline failure (#84 fallback). No call to reduce_linear_form, no de Weger lattice: the inequality above is the algorithm of this child.

Round record (unified #72/#82 field names, Dujella–Pethő semantics):

{"dimension": 2,                  # the Baker–Davenport approximation is 2-dimensional (kappa, mu); literal 2
 "scaling_C": q,                  # the convergent denominator (exact ZZ) — the round's scaling modulus
 "precision_bits": b,             # working ball precision that certified the round (ZZ)
 "b1_norm_lower": eps_lo,         # the certified exact rational lower bound for Lemma 5a's eps
 "new_bound": B_{r+1}}            # ZZ

In-memory values exact ZZ/QQ, tagged only in Evidence.to_dict() — the #72/#82 exactness rule.

Rounds and termination (exact rule, unchanged): from the current bound B_r, one application of the corollary produces B_{r+1}; iterate. If the new bound is >= the previous one, stop; the final bound is the previous one. The loop strictly decreases while it runs, so it terminates. Matveev and the Dujella–Pethő rounds run once per ray pair — and by the one-sided-orbit lemma there is exactly one ray per seed, so ray pairs are in canonical bijection with seed pairs. In this unit regime each equation has exactly one seed class (the norm-1 solutions are {±eps^k}, and any two norm--1 solutions differ by a norm-1 unit), hence there is exactly one ray pair; the per-ray-pair loop is real generality for #83, which inherits it unchanged. Illustrative chain on the classical fixture (verified against Sage 10.7; committed values come from the implementation's first verified run, #82): starting from M = 5e15: q = 46227103821174230, eps = 0.335, new bound 15; then M = 15: q = 518, eps = 0.105, new bound 3; then M = 3: q = 83, new bound 2; then a non-improving round — stop at final bound 2.

Finishing walk: convert the final per-pair index bounds to a certified stop height for equation 1's coordinates — H_stop := max over rays r1 of ceil(upper_ball(beta_{r1} * eps_1^(J(r1)))) with J(r1) := max over r2 of B_final(r1, r2) (each B_final >= j0, jz by construction, so every exceptional small index, the z^2 <= 2*T0 strip, and the general-regime x_1 = x_2 = 0 point are below the stop height) — then reuse the existing enumeration: walk _pell_solutions_stream(D1, N1), testing each shared value in equation 2 exactly as the current solver does, stopping when the stream's sup-norm level exceeds H_stop (beta*eps_1^J bounds both coordinates at index J). Collect matches, expand by the sign symmetries into full tuples (the walk's ± orbit structure already yields them, and streams the z = 0 seeds when N_1 is a square), and return the sorted list. This finishing step is today's loop with a certified stopping rule replacing the hardcoded 400/10^40 bounds.

Certified precision and rounding

Per #72: all comparisons on real quantities go through ball enclosures with one-sided conclusions only (a bound is used only in the direction its enclosure certifies — upper bounds where upper bounds are consumed, lower where lower; the roundings for A, eps_lo, W, H_stop, kappa_M, and the A_j are pinned above); precision doubling on inconclusive comparisons up to the cap; no floats anywhere in the decision path. The only #72 surface this child consumes: exact-value ball evaluation, the escalation/cap discipline, and the record conventions.

Failure and fallback semantics

• Regime detection fails (nonunit N, dependent fields, empty second seed set, non-concrete data, shared unknown not the B-coordinate of both): fall through to the current code path unchanged — exact fallback outputs are Simultaneous-Pell solver integration and partial fallback boundaries #84's contract; this child only asserts that its branch declines.
• First equation's seed set empty: the existing proved-empty branch answers (unchanged).
• Certified logarithmic-bound and LLL-reduction primitives #72 precision cap reached before a round certifies, or before the bound stabilizes: pipeline failure, fall through to the current partial path (Simultaneous-Pell solver integration and partial fallback boundaries #84).
• A cas-failure from Complete seed-orbit enumeration for generalized Pell equations #80/PARI propagates exactly as it does today.
• The proved branch never returns a partial result: it either succeeds with the contract above or declines entirely.

Acceptance criteria

1. Classical fixture, end-to-end, written in full. solve("x^2 - 3*z^2 = 1; y^2 - 8*z^2 = 1") returns exactly the 12 tuples, ordered (x, z, y) (first appearance), sorted ascending:

   (-2, -1, -3), (-2, -1, 3), (-2, 1, -3), (-2, 1, 3),
   (-1, 0, -1), (-1, 0, 1), (1, 0, -1), (1, 0, 1),
   (2, -1, -3), (2, -1, 3), (2, 1, -3), (2, 1, 3)
   

   — all sign combinations of (1, 0, 1) and (2, 1, 3); this is the set already returned by today's partial search and verified complete by brute force to z <= 10^6; this child upgrades the same list to kind="finite-complete", completeness="proved", scope="all integer solutions", S.complete is True.

2. Ray-decomposition unit tests (the new lemma, pinned). one_sided_rays(seed_orbits(13, 27)) returns exactly the four rays of the verified table above — seeds [(-12, 3), (12, 3), (-40, 11), (40, 11)] with (sign, shift, base) equal to (-1, 1, (768, 213)), (+1, 0, (12, 3)), (-1, 1, (220, 61)), (+1, 0, (40, 11)) in that order — and a test reproduces the exactly-one coverage check against an in-test brute-force enumeration of all signed solutions with |z| <= 10^5: each sign-normalized representative with z != 0 occurs in exactly one ray at exactly one index, ray z-sets are pairwise disjoint, and every ray element passes the exact equation check with x >= 0 < z. one_sided_rays(seed_orbits(3, 1)) returns the single ray (+1, 1, (2, 1)) and the coverage test additionally pins the z = 0 handling ((±1, 0) in no ray). Ray order equals seed order.

3. Ordering variant. solve("y^2 - 8*z^2 = 1; x^2 - 3*z^2 = 1") returns the same solution set ordered (y, z, x), i.e. exactly the sorted list of all sign combinations of (1, 0, 1) and (3, 1, 2) — pinning the first-appearance ordering convention.

4. Bennett cross-check assertion. For (N1, N2) = (1, 1) the test asserts that the number of returned tuples with all coordinates strictly positive is at most three (Bennett1998); for the fixture it is exactly one, (2, 1, 3).

5. A negative-Pell instance as a required assertion. For solve("x^2 - 2*z^2 = -1; y^2 - 3*z^2 = 1") (N1 = -1 with nonempty seed set — (1, 1) solves it, and its single ray has base (1, 1), shift 0; D1*D2 = 6 is not a square) assert completeness == "proved" and that the returned list equals the independent brute-force enumeration of all solutions with |z| <= 10^6 computed inside the test; the literal list is committed from the first verified run and asserted thereafter. (The eight tuples (±1, ±1, ±2), ordered (x, z, y), must be among the solutions: (±1)^2 - 2*(±1)^2 = -1, (±2)^2 - 3*(±1)^2 = 1.)

6. Firewall and serialization. Every returned point passes solvers._verified against the original system; as_dict() stays JSON-serializable with exact tagged values.

7. Evidence. The evidence tuple contains the two theorem records (reference keys BugeaudMignotteSiksek2006 and DujellaPetho1998, both resolving in data/references.bib — add the latter per the parent, coordinating with Rigorous de Weger completion of the exponential searches #55's children) and the computation record; all kinds from EVIDENCE_KINDS. Reduction-round entries carry exactly the five unified keys {"dimension", "scaling_C", "precision_bits", "b1_norm_lower", "new_bound"} with dimension == 2 and exact ZZ/QQ in-memory values; a test asserts the key set literally.

8. Regime boundary (light). solve("x^2 - 2*z^2 = 1; y^2 - 8*z^2 = 1") (dependent fields, D1*D2 = 16 square) still returns completeness == "partial" — the proved branch declines; the byte-exact fallback assertion is Simultaneous-Pell solver integration and partial fallback boundaries #84's fixture.

9. Doctest flip. _solve_simultaneous_pell's EXAMPLES assert the new proved behavior; make doctest is clean. New functions (including one_sided_rays and every private helper) carry Sage-convention docstrings with INPUT/OUTPUT and EXAMPLES passing sage -t.

10. Capability. CAPABILITIES["simultaneous-pell"].completeness == frozenset({"partial", "proved"}), .kinds == frozenset({"witness", "finite-complete", "empty"}), goals unchanged; the solver-contract suite's constraints (proved implies complete + evidence + scope) pass on the new results. make test, make doctest, make coverage (docstring coverage stays 100%), make registry-docs clean.

Out of scope

• Nonunit right-hand sides (Extension to supported nonunit right-hand sides #83, which reuses one_sided_rays and the per-ray-pair analysis verbatim).
• The final certificate schema and its exact JSON contract (Evidence and certificate integration for simultaneous Pell bounds #82) — this child's computation-record payload is provisional in structure (its round entries already use the unified field names).
• The fallback-boundary enumeration and byte-exact fallback fixtures (Simultaneous-Pell solver integration and partial fallback boundaries #84); any change to the partial path's strings, bounds, or plain-string evidence.
• Seed-enumeration internals (Complete seed-orbit enumeration for generalized Pell equations #80) — one_sided_rays is additive.

References

BugeaudMignotteSiksek2006 (exists; Theorem 9.4 transcribed above from arXiv:math/0403046 — cite as "Matveev's theorem in the form of [BugeaudMignotteSiksek2006, Thm 9.4]", the same convention as #73), DujellaPetho1998 (added per the parent; Lemma 5a stated in full above, transcribed from the author-hosted preprint and cross-checked against arXiv:1510.05579 — coordinate with #55's children so the bibliography gains the entry exactly once), Bennett1998 (exists; the at-most-three-positive-solutions cross-check), Anglin1996 (exists; family reference). The Nagell bounds live in #80.