diff --git a/.github/workflows/c8-window-coefficient-d6.yml b/.github/workflows/c8-window-coefficient-d6.yml new file mode 100644 index 00000000..ab9554af --- /dev/null +++ b/.github/workflows/c8-window-coefficient-d6.yml @@ -0,0 +1,79 @@ +name: C8 window coefficient d6 + +'on': + pull_request: + paths: + - 'frontiers/c8_window_coefficient_d6_20261001/**' + - '.github/workflows/c8-window-coefficient-d6.yml' + - 'imports/lifetime_parent_20260925/UNIFORM_MATRIX_CAP_AND_LIFETIME.md' + - 'coefficients/side24_v1/PROOF.md' + workflow_dispatch: + +permissions: + contents: read + +jobs: + certify-and-reject-mutants: + runs-on: ubuntu-latest + timeout-minutes: 90 + steps: + - uses: actions/checkout@11d5960a326750d5838078e36cf38b85af677262 + with: + persist-credentials: false + - uses: actions/setup-python@a26af69be951a213d495a4c3e4e4022e16d87065 + with: + python-version: '3.11.16' + - name: Verify sources and pins, replay both modes, reject mutants + env: + PYTHONDONTWRITEBYTECODE: '1' + run: | + set -euo pipefail + python -B -S - <<'PY' + import concurrent.futures, hashlib, json, pathlib, subprocess, sys + root = pathlib.Path('frontiers/c8_window_coefficient_d6_20261001').resolve() + manifest_doc = json.loads((root / 'SOURCE_FILES.json').read_text()) + manifest = manifest_doc['files'] + present = sorted(p.name for p in root.iterdir()) + expected_names = sorted([e['path'] for e in manifest] + ['SOURCE_FILES.json']) + if present != expected_names: + raise ValueError('packet tree differs from manifest: ' + repr(present)) + for entry in manifest: + path = root / entry['path'] + if pathlib.PurePosixPath(entry['path']).name != entry['path'] or path.is_symlink(): + raise ValueError('flat regular source required: ' + entry['path']) + data = path.read_bytes() + if len(data) != entry['bytes'] or hashlib.sha256(data).hexdigest() != entry['sha256']: + raise ValueError('source identity mismatch: ' + entry['path']) + # main-resident pins: the consumed sources must be the cited bytes + for pin in manifest_doc['pins_on_main']: + if not pin.get('workflow_checked'): + continue + p = pathlib.Path(pin['path']) + if p.is_symlink() or not p.is_file(): + raise ValueError('pinned source missing: ' + pin['path']) + data = p.read_bytes() + blob = hashlib.sha1(b'blob %d\0' % len(data) + data).hexdigest() + if blob != pin['git_blob']: + raise ValueError('pinned source drifted: ' + pin['path'] + ' ' + blob) + script = str(root / 'window_coefficient_d6.py') + expected = (root / 'RESULTS.json').read_bytes() + mutants = ['image-shells', 'd5-floor', 'sandwich-swap', 'window-scale', 'trace-slope', 'moment-recurrence', 'remainder-dropped', 'parts-sign', 'schur-sign', 'tau-cross', 'sphere', 'tilt-gamma', 'layer-dropped', 'cumulative-shift'] + def mutant_rc(args): + flags, m = args + return m, flags, subprocess.run([sys.executable, *flags, script, '--mutant', m], + capture_output=True, timeout=1500).returncode + for flags in (['-B', '-S'], ['-B', '-O', '-S']): + out = subprocess.run([sys.executable, *flags, script], + check=True, capture_output=True, timeout=1500).stdout + if out != expected: + raise RuntimeError('output differs from RESULTS.json under ' + ' '.join(flags)) + doc = json.loads(out) + if not (doc['passed'] and len(doc['checks']) == 19 and all(doc['checks'].values())): + raise RuntimeError('certification did not pass under ' + ' '.join(flags)) + with concurrent.futures.ThreadPoolExecutor(max_workers=4) as pool: + for m, fl, rc in pool.map(mutant_rc, [(flags, m) for m in mutants]): + if rc != 1: + raise RuntimeError('mutant ' + m + ' not rejected under ' + ' '.join(fl)) + print(json.dumps({'scientific_effect': 'NONE', 'mathematical_acceptance': False, 'certified_enclosures': True, 'C8': 'OPEN'})) + PY + git diff --exit-code diff --git a/frontiers/c8_window_coefficient_d6_20261001/NOTE.md b/frontiers/c8_window_coefficient_d6_20261001/NOTE.md new file mode 100644 index 00000000..654c0aca --- /dev/null +++ b/frontiers/c8_window_coefficient_d6_20261001/NOTE.md @@ -0,0 +1,379 @@ +# The compact-window coefficient `c_{B,K}` in `d = 6`: reference value and transfer to the torus for every `L ≥ 10` and every frame + +**Object** `CL-C8-WINDOW-COEFFICIENT-D6-20261001-v1` · **scientific effect NONE** · declarative record, `executed: false` · +base `main 124c37d` · claim: Math-#213 comment 5930520277. + +Companions (not on `main`): +- Math-#213: the `d = 5` record. Its `m = 4` layers become here one layer recursion for every `m ≤ 5`. +- Math-#209: the total-trace reduction and the Taylor quadrature. +- Math-#205: Lemma S. +- Math-#197: the interval toolkit, Owen's `T`, and the `d = 3` reference closed form. +- Math-#204: the certified `D_5`, `c_{6,ref}`, `Z_5`, the `d = 6` image-bound constant, and the cumulative layer for Owen-type terms. +- Math-#202 and #201: the closed forms of `D_4` and `D_3`. + +Catalog entry C8 stays OPEN; no register, GRAPH, STATUS or catalog surface is touched. Same GitHub account as every lane: zero +organizational-independence credit. Claude reads of this record count for nothing. The author will not merge. + +## 0. Statement + +Let `f` be the unit-variance stationary Gaussian field on `R⁶/(LZ⁶)` with the parent kernel `K_L` of [LP] §1. Let `B = [b_−, b_+]` +be a birth window and `K = [k_−, k_+]` a gap window, compact and of positive length with a positive gap floor: +`−∞ < b_− < b_+ < ∞` and `0 < k_− < k_+ < ∞`. This is [LP] Theorem B's hypothesis (lines 16 and 37 of the pinned [LP]). Let +`c^(6)_{B,K}` be the leading coefficient of [LP] Theorem B for the pairs with birth in `B` and scaled gap in `K` ([LP] (11.3), §15; +both stated for every `d`). Windows outside that hypothesis are treated in §7. + +**Proposition (d = 6 window coefficient).** +- **Reference kernel.** `c^(6)_{B,K} = c_{6,ref} F^(6)_B G_K`, with the birth factor `F^(6)_B` given by the one-dimensional integral + (2.3) below. +- **Torus.** For every `L ≥ 10` and every frame, `c^(6)_{B,K}` lies in the interval computed by `window_coefficient_d6.py` from the + entry box of radius `E^(6)_10 ≈ 4.610·10⁻⁹`. + +On the C8 band `B = [0, 1]`, `K = [1/2, 2]`: + +``` +c^(6)_{B,K}(reference) = 0.000004147268276047706822538464047034235965965259540700346…, +c^(6)_{B,K}(torus, every L ≥ 10, every frame) ∈ [0.0000041472119578916668252, 0.0000041473245949635375804] (relative half-width 1.36e-5), +c^(6)_{B,K}(torus, L = 24) = the reference digits to about 50 places (quadrature-limited). +``` + +The same holds on every window of §5, with relative half-width at most `1.6·10⁻⁵` at `L ≥ 10`. The full window returns +`c_{6,ref}`, Math-#204's value, and narrows its `D_5` enclosure from width `3.9·10⁻³⁷` to `8.8·10⁻⁵⁵`. + +## 1. The object and the reference law + +The setup is Math-#209 §1 with `|S⁵| = π³` and frame `(u, w_1, …, w_5)`: +- `G = ∇f ∈ R⁶`; +- `V = (f_uu, f_uw1, …, f_uw5)`; +- `A` is the `5×5` transverse Hessian; +- `t = f_uuu`. + +``` +c^(6)_{B,K} = 144 ∫_{S⁵} p_G(0) p_V(0) I_B(u) J_K(u) dσ(u), I_B(u) = E[1{f ∈ B} det(A)² 1{A ≺ 0} | V = 0]. (1.1) +``` + +For the reference kernel, `Cov V = diag(3, 1, 1, 1, 1, 1)`, `Cov G = I` and `τ² = 6`. For the 16-vector `(f, A_11, A_12, …, A_55) | V = 0`: + +``` +C'_ref: Var f = 2/3, Cov(f, A_ii) = −2/3, Var A_ii = 8/3, Cov(A_ii, A_jj) = 2/3, Var A_ij = 1 (i ≠ j), all else 0; +A | f = b, V = 0 ~ −b I + Q, Q the 5×5 GOE (diagonal variance 2, off-diagonal 1). (1.2) +``` + +This is checked by rule `REFERENCE_EXACT`. + +## 2. The reference birth integral + +Math-#209 (2.1) in general `m` puts the window on the total trace `T = −Σ μ_i`. For `m = 5`, use the ordered sector `μ_5 = −a`, +`μ_4 = −a − p_1`, …, `μ_1 = −a − p_1 − … − p_4` and the coordinates + +``` +T = 5a + 4p_1 + 3p_2 + 2p_3 + p_4, w′ = 4p_1 + 3p_2 + 2p_3 + p_4, w = 3p_2 + 2p_3 + p_4, s = 2p_3 + p_4, q = p_4 +(Jacobian 1/120, region T > w′ > w > s > q > 0). +``` + +The `5!` sectors cancel the Jacobian, and the exponent decouples exactly: +`|μ|²/4 − T²/32 = 3T²/160 + w′²/80 + w²/48 + s²/24 + q²/8`. The coefficients `1/8, 1/24, 1/48, 1/80` are `1/(8·(1, 3, 6, 10))`, as in +Math-#204. With Mehta's `Z_5 = 46080 π^{3/2}`: + +``` +I^ref_B = (√(3/8)/(46080 π^{3/2})) ∫_{T>w′>w>s>q>0} P e^{−3T²/160 − w′²/80 − w²/48 − s²/24 − q²/8} W_B(T) dq ds dw dw′ dT, +P = Π μ_i² · Π_{i 1`, the script evaluates +`T_Owen(h, a)` by the reflection `T(h, a) = ½[Φ(h) + Φ(ah)] − Φ(h)Φ(ah) − T(ah, 1/a)` and Math-#197's Owen series (1400 terms). + +Rule `INNER_LAYERS` requires five things: +- (2.2) holds coefficient by coefficient; +- at each layer, the set of kinds that the derivation generates and that come out zero is exactly the one stated above. A missing kind + fails the rule, so each cancellation is derived, not assumed; +- the outer exponent is `3/160`; +- the same recursion at `m = 4` reproduces Math-#213's stated `F_4` exactly; +- the interval value of (2.2) at `T = 3, 6, 10`, with (2.4) for the Owen-type kinds, agrees with a direct floating-point 4-D + Gauss–Legendre integral of (2.1) (20 nodes per layer) to `10⁻¹⁰` relative; the observed agreement is about `10⁻¹³`. + +**Controls.** +- **The cumulative layer against Owen's `T` (rule `OWEN_EXACT`).** The quadrature carries each `K_{α,g}` as a certified cumulative + layer (§3). At the cell centres `T = 1/4` and `81/4`, and at the cut `T = 100`, its value meets (2.4), both to `10⁻⁵⁰`; for + example `K_{1/80,1/16}(100) = 20.57651203962183075809061835089575319152259889923…`. +- **Full window (rule `D5_CERTIFIED`).** For `B = R`, (2.3) gives + `D_5 = 44.13065187507413271036202489540263314878223806307704201…` (width `8.8·10⁻⁵⁵`). It meets Math-#204's certified interval, + which has width `3.9·10⁻³⁷`. +- **Same code at `m = 2, 3, 4` (rule `CROSSCHECK_M2_M4`).** At `m = 2` it returns `29/6 − √6` and Math-#197's Owen-`T` closed form on + three windows. At `m = 3` it returns Math-#201's `D_3` and Math-#209's band factor `F^(4)_{[0,1]}`. At `m = 4` it returns + Math-#202's closed form of `D_4` and Math-#213's band factor `F^(5)_{[0,1]}`. +- **Independent outer quadrature (rule `FLOAT_CONTROL`).** The `T`-integral (2.3) by Gauss–Legendre (four panels of 48 nodes on + `(0, 80)`), with `F_5` evaluated in floating point and `K_{α,g}` by Simpson's rule, is independent of the Taylor quadrature and of + the cumulative layer. It agrees with (2.3) on `[0, 1]` and `[1, 2]` to `10⁻⁸` relative. + +## 3. Certified quadrature + +The quadrature is Math-#209 §3: + +| parameter | value | +|---|---| +| Taylor order | `48` | +| cell width | `1/2` | +| range | `(0, 100)` | +| Cauchy radius | `6` | + +The Owen-type kinds enter as a certified cumulative layer, the method of Math-#204 §3. At each cell centre `s_0`: +- the value of `K_{α,g}(s_0)` is the certified integral of `k = e^{−αw²} E_g(w)` over the earlier full cells, plus its integral over + the left half of the current cell, plus a one-sided Cauchy remainder; +- the higher Taylor coefficients of `K` at `s_0` are those of `k` divided by `n`; +- the disc bound of `K` is `|K(s_0)| + ρ M_k`. + +Beyond the cut, the integrand is bounded using `|K_{α,g}| ≤ π/(4√(αg))`. The Cauchy remainders sum to `6.2·10⁻⁵⁵` and the tail is +`5.4·10⁻⁶⁵`. Rule `QUAD_COARSE` runs the same quadrature at order `8` with cells of width `2`; its enclosure, remainder included, is +wide (`[7.95, 80.31]`) but must contain `D_5`. + +## 4. The torus in `d = 6` + +**Lemma S** (Math-#205 §2) applies verbatim to the 16-vector `(f, A)`, whose `det(A)²` has degree `10`. If +`(1 − ε)C'_ref ≤ C' ≤ (1 + ε)C'_ref`, then + +``` +(1 − ε)^13/(1 + ε)^8 · I^ref_{B/√(1−ε)} ≤ I'_B ≤ (1 + ε)^13/(1 − ε)^8 · I^ref_{B/√(1+ε)}. (4.1) +``` + +**Eigenvalue floor.** In general `m`, the eigenvalues of `C'_ref` are: +- `1` on the off-diagonal `A_ij`; +- `2` on the traceless diagonal; +- the two eigenvalues of `[[2/3, −(2/3)√m], [−(2/3)√m, (2m+6)/3]]` on `(f, tr A/√m)`. That block has determinant `4/3`. + +For `m = 5` this gives `λ_min = (6 − √(92/3))/2 = 0.23113`. The `d = 5` floor `1/4` therefore fails. Rule `LAMBDA_FLOOR` checks, by +Sylvester's criterion in rationals, that `C'_ref − (23/100) I` and `C'_ref − 0.231 I` are positive definite and that +`C'_ref − 0.2312 I` and `C'_ref − I/4` are not. Hence `ε := (100/23)‖C' − C'_ref‖_F`. + +**Image bound.** The shell `|n|_∞ = j` holds `(2j+1)⁶ − (2j−1)⁶ = 384j⁵ + 320j³ + 24j ≤ 728j⁵` points, with `|n|² ≤ 6j²`. So +`Σ_{n≠0}|n|⁶ e^{−L²|n|²/2} ≤ 157248 Σ_j j¹¹ e^{−L²j²/2} ≤ 314496 e^{−L²/2}` for `L ≥ 10` (successive terms have ratio at most +`2¹¹ e^{−3L²/2} < 1/2`); this is Math-#204's `E_6` constant. Every covariance entry of the 6-jet of `K_L`, in every frame, is within + +``` +E^(6)_L = 314496 (76 L⁶ + 15) e^{−L²/2}: E^(6)_10 ≈ 4.6100e-9, E^(6)_12 ≈ 3.8398e-18, E^(6)_24 ≈ 3.8272e-110 (rule IMAGE_BOUND) +``` + +of its reference value. The contraction bound, `|D^q φ(0)| ≤ 15` and the normalization are dimension-free (Math-#197 §3). + +**Interval evaluation.** `Space6(E)` widens by `±E` every entry of the 22-variable even block `(f, f_uu, f_uw1, …, f_uw5, A)` (253 +entries; all but `Var f = 1`, which stays exact) and the 28 entries of the odd block `(f_u, f_w1, …, f_w5, t)`. Over that box it +evaluates the following by interval linear algebra (adjugate inverses, Laplace determinants): +- `p_G(0) = (2π)^{−3} det(Cov G)^{−1/2}`; +- `p_V(0) = (2π)^{−3} det(Cov V)^{−1/2}`; +- `τ²`; +- the `16×16` Schur complement `C'`. + +Rule `EPSILON` checks the resulting values: + +``` +ε(L ≥ 10) ≤ 7.072e-7, ε(L ≥ 12) ≤ 5.890e-16, ε(L = 24) ≤ 3.4e-57 (arithmetic-limited), ε(E = 0) ≤ 1.1e-58; +τ²(L ≥ 10) ∈ [5.999999511, 6.000000489], p_G p_V (L ≥ 10) ∈ [0.000009383398295622, 0.000009383398785879] (reference (2π)^{−6}/√3). +``` + +Then `c^(6)_{B,K} ∈ 144 · π³ · [p_G p_V] · (4.1)(I^ref, ε) · J_K(τ²)`. + +**Exact tests of (4.1).** +- **Rule `SANDWICH_EXACT_SCALING`.** For `C' = (1 + η)C'_ref` with `η = ±1/50`, on `[0, 1]`, `[3, 7/2]` and `[−7/2, −3]`, the exact + value `(1 + η)⁵ I^ref_{B/√(1+η)}` lies inside its own sandwich. +- **Rule `TILTED_EXACT`.** This is Math-#209's tilted pair `(v, μ, σ²) = (17/25, 33/34, 133/5100)`, which is non-scalar, with + `ε = (100/23)·√14/150 = 0.1085`. Its exact value is evaluated through the same reduction. On `B = R` it equals `D_5`. On `[0, 1]` + and `[1, 2]` it exceeds the reference by `29.9 %` and `7.0 %`, and lies inside (4.1). + +## 5. Values (`RESULTS.json`) + +`F^(6)_B = I^ref_B/D_5`, `G_K` is Math-#197's, and `c^(6,ref) = c_{6,ref} F^(6)_B G_K`, with +`c_{6,ref} = Γ(7/6) (3/2)^{1/3} D_5/(64√3 π^{7/2})` (Math-#204). Displayed intervals are rounded outward (lower endpoints down, upper +endpoints up) and contain the certified intervals of `RESULTS.json`. Every `…`-terminated digit string is a common prefix of a +certified interval's endpoints. Both properties were checked mechanically before publication. + +| `B` | `K` | `F^(6)_B` | `c^(6)` reference | torus, every `L ≥ 10`, every frame | +|---|---|---|---|---| +| `[0,1]` | `[1/2,2]` | `0.00800194436361…` | `0.0000041472682760477068225384…` | `[0.0000041472119578916668252, 0.0000041473245949635375804]` | +| `[0,1]` | `(0,∞)` | `0.00800194436361…` | `0.0000615579867872192435333128…` | `[0.000061557165243731316444, 0.000061558808341600505797]` | +| `[−1,1]` | `[1/2,2]` | `0.00802638034698…` | `0.0000041599330202636062388615…` | `[0.0000041598765100768776358, 0.0000041599895312130049122]` | +| `[−2,2]` | `[1/2,2]` | `0.21153125452829…` | `0.0001096329618694244702024245…` | `[0.00010963146730285680056, 0.00010963445645615433875]` | +| `[0,∞)` | `[1/2,2]` | `0.99997555972228…` | `0.0005182699013148570581656579…` | `[0.00051826204265890159945, 0.00051827776008865107624]` | +| `(−∞,0]` | `[1/2,2]` | `0.00002444027771…` | `1.2666969903564286844701424489…·10⁻⁸` | `[1.2666777831135967712e-8, 1.2667161978872683484e-8]` | +| `R` | `[1/2,2]` | `1` | `0.0005182825682847606224525026…` | `[0.00051827470943673273541, 0.00051829042725062994893]` | +| `R` | `(0,∞)` | `1` | `0.0076928786292368482786666312…` (= `c_{6,ref}`) | `[0.0076927637782619387375, 0.0076929934819071653012]` | + +`RESULTS.json` lists all twelve window pairs, and each `L ≥ 10` enclosure contains its reference value. At `L = 24` every enclosure +agrees with the reference to about 50 digits. The full window carries Math-#204's 40 certified digits of `c_{6,ref}` +(rule `C6_CONSISTENT`). + +**How the band's share moves with dimension.** + +| `d` | height-window share `F^(d)_{[0,1]}` (every `K`) | C8 band's share of `c_{d,ref}`, `F^(d)_{[0,1]} G_{[1/2,2]}` | +|---|---|---| +| 2 | `48.3 %` | `3.25 %` | +| 3 | `34.6 %` | `2.33 %` | +| 4 | `15.4 %` | `1.04 %` | +| 5 | `4.37 %` | `0.295 %` | +| 6 | **`0.800 %`** | **`0.0539 %`** | + +The first column is the share of birth heights in `B = [0, 1]`; by the factorization `c_{B,K} = c_{d,ref} F_B G_K` it does +not depend on the gap window. The second is the C8 band's share of the whole coefficient, with the dimension-free gap factor +`G_{[1/2,2]} = 6.74 %`. Both are ratios of certified values from the separate packets, not a theorem about birth mass. The negative half-line carries `0.0024 %` of the birth heights in `d = 6`. + +## 6. Rules, mutants, verification + +`python3 -B -S window_coefficient_d6.py` (also `-B -O -S`; about two minutes) prints `RESULTS.json`. It exits `0` only if all +nineteen rules hold: + +| rule | what it requires | +|---|---| +| `REFERENCE_EXACT` | the exact `C'_ref` and conditional law, including the tilted covariance at `(2/3, 1, 0)` | +| `LAMBDA_FLOOR` | the floor `23/100` certified; `0.2312` and `1/4` refused | +| `INNER_LAYERS` | (2.2) exact, the per-layer cancellations, `m = 4` against Math-#213, and the float 4-D integral | +| `OWEN_EXACT` | the cumulative layer meets Owen's `T` (2.4) at three points for each kind | +| `D5_CERTIFIED` | the full window meets Math-#204's `D_5`, width below `10⁻⁵⁰` | +| `QUAD_COARSE` | the coarse quadrature with its remainder contains `D_5` | +| `CROSSCHECK_M2_M4` | the `m = 2, 3, 4` runs reproduce the `d = 3, 4, 5` results | +| `FLOAT_CONTROL` | the independent outer float quadrature agrees with (2.3) | +| `IMAGE_BOUND` | the `d = 6` image-bound constant and values | +| `EPSILON` | the certified `ε` values | +| `SANDWICH_EXACT_SCALING` | exact scalar perturbations lie inside (4.1) | +| `TILTED_EXACT` | the tilted non-scalar perturbation lies inside (4.1) | +| `FACTORIZATION` | the pipeline at `E = 0` reproduces `c_{6,ref} F^(6)_B G_K` on the band and the full window | +| `C6_CONSISTENT` | agreement with Math-#204's `c_{6,ref}` | +| `NESTING` | `L ≥ 10` enclosures contain the reference; the band meets the `10³⁰`-times narrower `L = 24` enclosure | +| `WIDTHS` | relative half-width below `2·10⁻⁵` at `L ≥ 10` | +| `SANDWICH_ORDERED` | every bracket has lower ≤ upper | +| `LIBRARY_EXACT` | exp, sqrt and negation against exact rational brackets | +| `PINNED` | the pinned digits | + +Fourteen mutants exit `1` in both interpreter modes: + +| mutant | change | rejected by | +|---|---|---| +| `image-shells` | `384j⁵` shells (the leading term only) | `IMAGE_BOUND` | +| `d5-floor` | the `d = 5` floor `1/4` | `LAMBDA_FLOOR`, `TILTED_EXACT` | +| `sandwich-swap` | the two rescalings exchanged | `SANDWICH_EXACT_SCALING`, `SANDWICH_ORDERED` | +| `window-scale` | `B·s` for `B/s` | `SANDWICH_EXACT_SCALING`, `SANDWICH_ORDERED` | +| `trace-slope` | the window factor's centre at `T/(m+2)` instead of `T/(m+3)` | `CROSSCHECK_M2_M4`, `FLOAT_CONTROL`, `PINNED` | +| `moment-recurrence` | `j` for `j + 1` in the moment recurrence | the recursion itself (see below) | +| `remainder-dropped` | the Cauchy remainders omitted, in the quadrature and in the cumulative layer | `QUAD_COARSE` | +| `parts-sign` | the sign of the boundary term in the integration by parts | `INNER_LAYERS`, `D5_CERTIFIED`, `CROSSCHECK_M2_M4`, `QUAD_COARSE`, `FLOAT_CONTROL` and four more | +| `schur-sign` | the sign of the Schur complement | `REFERENCE_EXACT`, `LAMBDA_FLOOR`, `TILTED_EXACT` | +| `tau-cross` | `Cov(t, f_u)` dropped | `REFERENCE_EXACT`, `FACTORIZATION`, `C6_CONSISTENT`, `NESTING` | +| `sphere` | `8π²/3` (the `d = 5` sphere) for `π³` | `FACTORIZATION`, `C6_CONSISTENT`, `NESTING` | +| `tilt-gamma` | the `σ²` term dropped from the tilted exponent | `TILTED_EXACT` | +| `layer-dropped` | the zero `E_{3/16}` layer silently left out of the derivation | `INNER_LAYERS` | +| `cumulative-shift` | the cumulative layer's value at each centre without the left half cell | `OWEN_EXACT`, `D5_CERTIFIED`, `QUAD_COARSE`, `FLOAT_CONTROL` and four more | + +`moment-recurrence` breaks the `w`-layer cancellation, so an Owen-type kind would enter the `w′`-layer; the recursion refuses that +case and the script exits `1` through the `ValueError`. `layer-dropped` changes no value, since the omitted layer is zero, and fails +`INNER_LAYERS` alone. `remainder-dropped` is the only mutant that `QUAD_COARSE` alone catches: without remainders the coarse +quadrature gives about `44.130659`, which misses `D_5 = 44.130651875…`. + +The hosted workflow replays the manifest, the pins, both modes byte for byte, and the mutants. + +## 7. What this does not do + +- It encloses `c^(6)_{B,K}` for the reference kernel and the torus. It does not revalidate [LP] Theorem B and §15, which it consumes + at their scope. It does not review Math-#197, #199, #201, #202, #204, #205, #209 or #213. +- It does not go to `d ≥ 7`. At `m = 6` one more layer (`c = 1/120`) follows the `w′`-layer, so the three Owen-type kinds of (2.2) + would have to be integrated against a further Gaussian layer unless they cancel there; the recursion as written refuses that case. + The floor also falls further (`0.2064` at `m = 6`). +- **Windows outside [LP] Theorem B's hypothesis.** Theorem B assumes compact windows of positive length with a positive gap floor: + `−∞ < b_− < b_+ < ∞` and `0 < k_− < k_+ < ∞`. Three kinds of window fall outside it: + - a gap window that touches `k = 0` or `∞`, such as the tabulated `K = (0,∞)`; + - an unbounded `B`: `[0,∞)`, `(−∞,0]`, `R`; + - a singleton window, `b_− = b_+` or `k_− = k_+`. There (1.1) is zero and does not inherit Theorem B's positive coefficient. No + singleton window is tabulated. + + For these windows the values are the coefficient integrals (1.1) on those sets only. Theorem B's asymptotic statement is not + extended to them. +- It does not address `C`, `r_*` or `z_*`. C8 stays OPEN; no register surface is touched; `executed: false`. +- The replay workflow follows the repository's per-packet convention. On the workflow question raised on Math-#202/#204, the + xAI lane ruled keep (thread reply 4155307498, a process ruling only); this workflow is the same construction. + +## 8. Provenance + +**Pins on `main 124c37d` (workflow-checked; the same blobs as on `main 3e0a91b`, the base of the earlier packets):** +- [LP] `imports/lifetime_parent_20260925/UNIFORM_MATRIX_CAP_AND_LIFETIME.md`, blob `dfed3b8d`. +- SIDE24 `coefficients/side24_v1/PROOF.md`, blob `44b66f04`. + +**Companions, not on `main`:** +- **Math-#197** at `ebc93cd`: `window_coefficient.py` blob `43b4f51b`. The toolkit, including Owen's `T`, is reused verbatim through + Math-#205's copy. +- **Math-#205** at `aebb4d4`: `window_coefficient_d3.py` blob `ba24b380`. Its toolkit section and `I3_ref` are reused verbatim; + Lemma S is its §2. +- **Math-#199** at `617236d`: `cone_moment_d4.py` blob `4787d6bf`. The exact polynomial algebra is reused verbatim. +- **Math-#209** at `84a494a`: `window_coefficient_d4.py` blob `f77f7d2d`. The reduction and the quadrature; its band factor is a rule. +- **Math-#213** at `b22ff4b`: `window_coefficient_d5.py` blob `da68e6ec`. The `m = 4` layers, generalized here; its stated `F_4` + and band factor are rules; the tilted pair and the exact tests. +- **Math-#204** at `5915af4`: `cone_moment_d6.py` blob `c5e56694` (the cumulative layer, `Z_5`, the `E_6` constant) and + `RESULTS.json` blob `031c8bbe` (`D_5`, `c_{6,ref}`). +- **Math-#202** at `4512309`: `RESULTS.json` blob `52ce9b57` (the closed form of `D_4`, `Z_4 = 1536π`). +- **Math-#201** at `3abcc10`: `NOTE.md` blob `cbf35dea` (the closed form of `D_3`). + +No external numerical library is used. + +## 9. Current disposition + +- **Head / owner:** see the PR's disposition block. Branch `claude/c8-window-d6-20261001` (base `main 124c37d`); author lane + Anthropic / Claude (`session_017Mi3hxjaxV45x6zo6o1ee3`). +- **Claim:** + - the `m = 5` total-trace reduction and its closed inner layers (2.2), including the recurring `w`-layer cancellation and the three + surviving Owen-type kinds with their Owen-`T` form (2.4); + - the certified reference values of `c^(6)_{B,K}` and the sharpened `D_5`; + - Lemma S in `d = 6` with the floor `23/100`; + - the `d = 6` image bound and the certified `ε`; + - the torus enclosures for every `L ≥ 10` and every frame, and for `L = 24`, on the windows of §5. + + No theorem of [LP] or SIDE24 is proved or reviewed. +- **Completed review scopes:** see the PR's disposition block, which is the live record (none at opening). Requested nonauthor reads: + - Slice A: §§1–3; + - Slice B: §4; + - Slice C: §§5–7. +- **Unresolved finding IDs:** none. +- **Validation:** 19/19 rules in both modes, byte-identical output; 14/14 mutants rejected in both modes; workflow replayed + locally; hosted runs: see the PR's disposition block (pending at opening). +- **Next action:** nonauthor reads; amendments on this branch, recorded in `SOURCE_FILES.json`. The author will not merge. + +## 10. Revisions + +- **v1.1 (band-share wording; NOTE only).** Two caption checks on this PR, the xAI lane (Harper, comment 5932429158) and a + scoped check (comment 5932380111; its lane is not stated), found that the table of §5, headed as the band's share of the + birth mass, shows `F^(d)_{[0,1]}`, the share of birth heights in `B = [0, 1]` (the same for every `K`), while it read as the + C8 band's share of the whole coefficient. §5 now labels that column as the height-window share and adds the C8 band's share + of `c_{d,ref}`, `F^(d)_{[0,1]} G_{[1/2,2]}`. The script, `RESULTS.json`, the rules and the mutants are byte-unchanged; no + certified value changes. +- **v1.2 (nonauthor read W3, Grok Bot agent 12, Math-#222 comment 5974639810, on `afd1c29`; NOTE §9 only).** §9 said + "Completed review scopes: none yet" and "hosted run pending at opening", which contradicted §10 and the PR body. It now points + to the PR's disposition block as the live record. The read found no mathematical defect. The script, `RESULTS.json`, the + rules, the mutants and the workflow are byte-unchanged. diff --git a/frontiers/c8_window_coefficient_d6_20261001/RESULTS.json b/frontiers/c8_window_coefficient_d6_20261001/RESULTS.json new file mode 100644 index 00000000..df361c43 --- /dev/null +++ b/frontiers/c8_window_coefficient_d6_20261001/RESULTS.json @@ -0,0 +1,1059 @@ +{ + "object": "CL-C8-WINDOW-COEFFICIENT-D6-20261001-v1", + "scientific_effect": "NONE", + "certified": true, + "mutant": null, + "passed": true, + "checks": { + "REFERENCE_EXACT": true, + "LAMBDA_FLOOR": true, + "INNER_LAYERS": true, + "OWEN_EXACT": true, + "D5_CERTIFIED": true, + "QUAD_COARSE": true, + "CROSSCHECK_M2_M4": true, + "FLOAT_CONTROL": true, + "IMAGE_BOUND": true, + "EPSILON": true, + "SANDWICH_EXACT_SCALING": true, + "TILTED_EXACT": true, + "FACTORIZATION": true, + "C6_CONSISTENT": true, + "NESTING": true, + "WIDTHS": true, + "SANDWICH_ORDERED": true, + "LIBRARY_EXACT": true, + "PINNED": true + }, + "c8_band": { + "B": "[0,1]", + "K": "[1/2,2]", + "c6_reference": [ + "0.00000414726827604770682253846404703423596596525954070034639824564", + "0.00000414726827604770682253846404703423596596525954070034643496455" + ], + "c6_reference_digits": "0.000004147268276047706822538464047034235965965259540700346", + "F6_B": [ + "0.00800194436361801032969394391850483826245390449805141765824915", + "0.00800194436361801032969394391850483826245390449805141766768055" + ], + "c6_torus_every_L_ge_10_every_frame": [ + "0.00000414721195789166682527996969593734250810229025973108427195893", + "0.00000414732459496353758032632074194675180160024735277578702022620" + ], + "c6_torus_L_24": [ + "0.00000414726827604770682253846404703423596596525954070034641389811", + "0.00000414726827604770682253846404703423596596525954070034641930870" + ], + "relative_half_width_L_ge_10": "0.0000135798547330598413438969465331434904143707892200054166738590" + }, + "full_window": { + "c_6_ref": [ + "0.00769287862923684827866663128122739008115923860384547469301787", + "0.00769287862923684827866663128122739008115923860384547469515793" + ], + "c6_torus_every_L_ge_10": [ + "0.00769276377826193873754843295895462363077900905449282168623764", + "0.00769299348190716530111789455991061068027529783142325169698437" + ], + "c6_torus_L_24": [ + "0.00769287862923684827866663128122739008115923860384547469247376", + "0.00769287862923684827866663128122739008115923860384547469570134" + ], + "Math204_c_6_ref_digits": "0.0076928786292368482786666312812273900811" + }, + "windows": { + "B=[0,1] K=[1/2,2]": { + "F6_B": [ + "0.00800194436361801032969394391850483826245390449805141765824915", + "0.00800194436361801032969394391850483826245390449805141766768055" + ], + "G_K": [ + "0.067371733425641667348844418474402773431540422159891980011888", + "0.067371733425641667348844418474402773431540422159891980510232" + ], + "c6_ref": [ + "0.00000414726827604770682253846404703423596596525954070034639824564", + "0.00000414726827604770682253846404703423596596525954070034643496455" + ], + "c6_torus_L_ge_10": [ + "0.00000414721195789166682527996969593734250810229025973108427195893", + "0.00000414732459496353758032632074194675180160024735277578702022620" + ], + "c6_torus_L_24": [ + "0.00000414726827604770682253846404703423596596525954070034641389811", + "0.00000414726827604770682253846404703423596596525954070034641930870" + ] + }, + "B=[0,1] K=[0,inf]": { + "F6_B": [ + "0.00800194436361801032969394391850483826245390449805141765824915", + "0.00800194436361801032969394391850483826245390449805141766768055" + ], + "G_K": [ + "0.999999999999999999999999999999999999999999999999999999742690", + "1.00000000000000000000000000000000000000000000000000000025731" + ], + "c6_ref": [ + "0.0000615579867872192435333128633011935775854894957305721223585503", + "0.0000615579867872192435333128633011935775854894957305721224799086" + ], + "c6_torus_L_ge_10": [ + "0.0000615571652437313164446963530097601960285483230793002875187250", + "0.0000615588083416005057968496129960795417953004215464123351318739" + ], + "c6_torus_L_24": [ + "0.0000615579867872192435333128633011935775854894957305721223717085", + "0.0000615579867872192435333128633011935775854894957305721224667252" + ] + }, + "B=[-1,1] K=[1/2,2]": { + "F6_B": [ + "0.00802638034698093333367660145536158010017940095668308069741289", + "0.00802638034698093333367660145536158010017940095668308070608238" + ], + "G_K": [ + "0.067371733425641667348844418474402773431540422159891980011888", + "0.067371733425641667348844418474402773431540422159891980510232" + ], + "c6_ref": [ + "0.00000415993302026360623886153497724923290474031172667552606581315", + "0.00000415993302026360623886153497724923290474031172667552610223439" + ], + "c6_torus_L_ge_10": [ + "0.00000415987651007687763581696579077369314035702698871425185661214", + "0.00000415998953121300491218558738599084910371193212124106830027303" + ], + "c6_torus_L_24": [ + "0.00000415993302026360623886153497724923290474031172667552608151364", + "0.00000415993302026360623886153497724923290474031172667552608653169" + ] + }, + "B=[-1,1] K=[0,inf]": { + "F6_B": [ + "0.00802638034698093333367660145536158010017940095668308069741289", + "0.00802638034698093333367660145536158010017940095668308070608238" + ], + "G_K": [ + "0.999999999999999999999999999999999999999999999999999999742690", + "1.00000000000000000000000000000000000000000000000000000025731" + ], + "c6_ref": [ + "0.0000617459698414162610822205661065417428332558222042764622594893", + "0.0000617459698414162610822205661065417428332558222042764623751355" + ], + "c6_torus_L_ge_10": [ + "0.0000617451454915504333130231724311003198955926488804660648765076", + "0.0000617467942022168491463245979265059981391159548677189246749750" + ], + "c6_torus_L_24": [ + "0.0000617459698414162610822205661065417428332558222042764622726909", + "0.0000617459698414162610822205661065417428332558222042764623619261" + ] + }, + "B=[-2,2] K=[1/2,2]": { + "F6_B": [ + "0.211531254528299509573293674591419351905112071481084366205592", + "0.211531254528299509573293674591419351905112071481084366230625" + ], + "G_K": [ + "0.067371733425641667348844418474402773431540422159891980011888", + "0.067371733425641667348844418474402773431540422159891980510232" + ], + "c6_ref": [ + "0.000109632961869424470202424576454217383529718997491924025057349", + "0.000109632961869424470202424576454217383529718997491924025911771" + ], + "c6_torus_L_ge_10": [ + "0.000109631467302856800564348024547241581179482062165667716529414", + "0.000109634456456154338742626534227726729027595029008866908395379" + ], + "c6_torus_L_24": [ + "0.000109632961869424470202424576454217383529718997491924025471118", + "0.000109632961869424470202424576454217383529718997491924025497957" + ] + }, + "B=[-2,2] K=[0,inf]": { + "F6_B": [ + "0.211531254528299509573293674591419351905112071481084366205592", + "0.211531254528299509573293674591419351905112071481084366230625" + ], + "G_K": [ + "0.999999999999999999999999999999999999999999999999999999742690", + "1.00000000000000000000000000000000000000000000000000000025731" + ], + "c6_ref": [ + "0.00162728426737641558642679261074172682261655922721833083026445", + "0.00162728426737641558642679261074172682261655922721833083174717" + ], + "c6_torus_L_ge_10": [ + "0.00162726246384221324672680310426002465830265924136651345353101", + "0.00162730607119968786820143362932637674507404497275232550089304" + ], + "c6_torus_L_24": [ + "0.00162728426737641558642679261074172682261655922721833083061223", + "0.00162728426737641558642679261074172682261655922721833083139939" + ] + }, + "B=[0,inf] K=[1/2,2]": { + "F6_B": [ + "0.999975559722285306230577770084739863440655934966762142287783", + "0.999975559722285306230577770084739863440655934966762142329072" + ], + "G_K": [ + "0.067371733425641667348844418474402773431540422159891980011888", + "0.067371733425641667348844418474402773431540422159891980510232" + ], + "c6_ref": [ + "0.000518269901314857058165657988231844700527927034400388612796138", + "0.000518269901314857058165657988231844700527927034400388616795321" + ], + "c6_torus_L_ge_10": [ + "0.000518262042658901599452227191534953538091592014569386347170021", + "0.000518277760088651076237752568089954999806954616279782540477243" + ], + "c6_torus_L_24": [ + "0.000518269901314857058165657988231844700527927034400388614748454", + "0.000518269901314857058165657988231844700527927034400388614842747" + ] + }, + "B=[0,inf] K=[0,inf]": { + "F6_B": [ + "0.999975559722285306230577770084739863440655934966762142287783", + "0.999975559722285306230577770084739863440655934966762142329072" + ], + "G_K": [ + "0.999999999999999999999999999999999999999999999999999999742690", + "1.00000000000000000000000000000000000000000000000000000025731" + ], + "c6_ref": [ + "0.00769269061314672429737085669538652941189912384812206112818662", + "0.00769269061314672429737085669538652941189912384812206113460309" + ], + "c6_torus_L_ge_10": [ + "0.00769257576497880447884516942468839347243354897459753051485012", + "0.00769280546301001016099992557289478525858024402723587856794652" + ], + "c6_torus_L_24": [ + "0.00769269061314672429737085669538652941189912384812206112977574", + "0.00769269061314672429737085669538652941189912384812206113301325" + ] + }, + "B=[-inf,0] K=[1/2,2]": { + "F6_B": [ + "0.0000244402777146937694222299152601365593440650332378576909248771", + "0.0000244402777146937694222299152601365593440650332378576922223494" + ], + "G_K": [ + "0.067371733425641667348844418474402773431540422159891980011888", + "0.067371733425641667348844418474402773431540422159891980510232" + ], + "c6_ref": [ + "1.26669699035642868447014244894914307042332908355743099975116E-8", + "1.26669699035642868447014244894914307042332908355743107671895E-8" + ], + "c6_torus_L_ge_10": [ + "1.26667778311359677121805019456876993139551782190539604094940E-8", + "1.26671619788726834833416058122805700139734361318107440325190E-8" + ], + "c6_torus_L_24": [ + "1.26669699035642868447014244894914307042332908355743100452279E-8", + "1.26669699035642868447014244894914307042332908355743107194668E-8" + ] + }, + "B=[-inf,0] K=[0,inf]": { + "F6_B": [ + "0.0000244402777146937694222299152601365593440650332378576909248771", + "0.0000244402777146937694222299152601365593440650332378576922223494" + ], + "G_K": [ + "0.999999999999999999999999999999999999999999999999999999742690", + "1.00000000000000000000000000000000000000000000000000000025731" + ], + "c6_ref": [ + "1.88016090123981295774585840860669260114755723413562687989089E-7", + "1.88016090123981295774585840860669260114755723413562698119448E-7" + ], + "c6_torus_L_ge_10": [ + "1.88013283134258703263534266230158345460079895291171377539570E-7", + "1.88018897155140117968987015825421695053804187373129047872966E-7" + ], + "c6_torus_L_24": [ + "1.88016090123981295774585840860669260114755723413562688027928E-7", + "1.88016090123981295774585840860669260114755723413562698080591E-7" + ] + }, + "B=[-inf,inf] K=[1/2,2]": { + "F6_B": [ + "0.999999999999999999999999999999999999999999999999999999980007", + "1.00000000000000000000000000000000000000000000000000000002000" + ], + "G_K": [ + "0.067371733425641667348844418474402773431540422159891980011888", + "0.067371733425641667348844418474402773431540422159891980510232" + ], + "c6_ref": [ + "0.000518282568284760622452502689656334191958631267691224187106809", + "0.000518282568284760622452502689656334191958631267691224191105417" + ], + "c6_torus_L_ge_10": [ + "0.000518274709436732735419939372036899225790905969747605401131102", + "0.000518290427250629948921235909695767280376968589715914351220600" + ], + "c6_torus_L_24": [ + "0.000518282568284760622452502689656334191958631267691224189059173", + "0.000518282568284760622452502689656334191958631267691224189152793" + ] + }, + "B=[-inf,inf] K=[0,inf]": { + "F6_B": [ + "0.999999999999999999999999999999999999999999999999999999980007", + "1.00000000000000000000000000000000000000000000000000000002000" + ], + "G_K": [ + "0.999999999999999999999999999999999999999999999999999999742690", + "1.00000000000000000000000000000000000000000000000000000025731" + ], + "c6_ref": [ + "0.00769287862923684827866663128122739008115923860384547469088460", + "0.00769287862923684827866663128122739008115923860384547469729125" + ], + "c6_torus_L_ge_10": [ + "0.00769276377826193873754843295895462363077900905449282168623764", + "0.00769299348190716530111789455991061068027529783142325169698437" + ], + "c6_torus_L_24": [ + "0.00769287862923684827866663128122739008115923860384547469247376", + "0.00769287862923684827866663128122739008115923860384547469570134" + ] + } + }, + "worst_relative_half_width_L_ge_10": "0.0000151635933714533483754430081497209514534878920287109230210644", + "transfer": { + "E6_10": [ + "4.61004016248414645975787970457461267471096818105826323160455E-9", + "4.61004016248414645975787970457461267471096818105826323160465E-9" + ], + "E6_12": [ + "3.83984351875935630376295613878951364604532248646439537253821E-18", + "3.83984351875935630376295613878951364604532248646439537253825E-18" + ], + "E6_24": [ + "3.82723721533215540079690742342942909477794757569364110281421E-110", + "3.82723721533215540079690742342942909477794757569364110281424E-110" + ], + "eps_L_ge_10": "7.07105202135832893622218649146813585623969607352992740183931E-7", + "eps_L_ge_12": "5.88969539292884405715595934327997441785915753506542895121270E-16", + "eps_L_24": "3.30389298929128313516258830463237243054623666474150446589803E-57", + "eps_E0": "1.0018885758089488924916983071823529312690153579497674688728E-58", + "lambda_floor": "23/100", + "lambda_min_exact": "(6 - sqrt(92/3))/2", + "tau2_L_ge_10": [ + "5.99999951133572415951581316086662414082700633849481172810193", + "6.00000048866424965744021031419012523577572718983858024716357" + ], + "pG_pV_L_ge_10": [ + "0.00000938339829562281227262380077247003879226230957839471447386805", + "0.00000938339878587838443632888368325784368016030300906761621913107" + ] + }, + "reference": { + "C_prime": [ + [ + "2/3", + "-2/3", + "0", + "0", + "0", + "0", + "-2/3", + "0", + "0", + "0", + "-2/3", + "0", + "0", + "-2/3", + "0", + "-2/3" + ], + [ + "-2/3", + "8/3", + "0", + "0", + "0", + "0", + "2/3", + "0", + "0", + "0", + "2/3", + "0", + "0", + "2/3", + "0", + "2/3" + ], + [ + "0", + "0", + "1", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0" + ], + [ + "0", + "0", + "0", + "1", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0" + ], + [ + "0", + "0", + "0", + "0", + "1", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0" + ], + [ + "0", + "0", + "0", + "0", + "0", + "1", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0" + ], + [ + "-2/3", + "2/3", + "0", + "0", + "0", + "0", + "8/3", + "0", + "0", + "0", + "2/3", + "0", + "0", + "2/3", + "0", + "2/3" + ], + [ + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "1", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0" + ], + [ + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "1", + "0", + "0", + "0", + "0", + "0", + "0", + "0" + ], + [ + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "1", + "0", + "0", + "0", + "0", + "0", + "0" + ], + [ + "-2/3", + "2/3", + "0", + "0", + "0", + "0", + "2/3", + "0", + "0", + "0", + "8/3", + "0", + "0", + "2/3", + "0", + "2/3" + ], + [ + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "1", + "0", + "0", + "0", + "0" + ], + [ + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "1", + "0", + "0", + "0" + ], + [ + "-2/3", + "2/3", + "0", + "0", + "0", + "0", + "2/3", + "0", + "0", + "0", + "2/3", + "0", + "0", + "8/3", + "0", + "2/3" + ], + [ + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "0", + "1", + "0" + ], + [ + "-2/3", + "2/3", + "0", + "0", + "0", + "0", + "2/3", + "0", + "0", + "0", + "2/3", + "0", + "0", + "2/3", + "0", + "8/3" + ] + ], + "A_given_f_mean": [ + "-1", + "0", + "0", + "0", + "0", + "-1", + "0", + "0", + "0", + "-1", + "0", + "0", + "-1", + "0", + "-1" + ], + "tau2": "6", + "det_Cov_V": "3", + "D_5_quadrature": [ + "44.1306518750741327103620248954026331487822380630770420119935", + "44.1306518750741327103620248954026331487822380630770420128758" + ], + "D_5_Math204": [ + "44.1306518750741327103620248954026331485852349159356435082422", + "44.1306518750741327103620248954026331489792412102184405166029" + ], + "coarse_quadrature": [ + "7.95875990612490307073900695599676371983111103705277410492199", + "80.3025585999501934825390870028869793108150551062709929593314" + ], + "F_5_pieces": { + "('K', Fraction(1, 30), Fraction(1, 24))": [ + "353396736/625", + "0", + "27240192/625", + "0", + "-3677184/3125", + "0", + "1085952/78125", + "0", + "-21504/390625", + "0", + "768/9765625" + ], + "('K', Fraction(1, 30), Fraction(1, 6))": [ + "706793472/625", + "0", + "54480384/625", + "0", + "-7354368/3125", + "0", + "2171904/78125", + "0", + "-43008/390625", + "0", + "1536/9765625" + ], + "('K', Fraction(1, 80), Fraction(1, 16))": [ + "265047552/625", + "0", + "20430144/625", + "0", + "-2757888/3125", + "0", + "814464/78125", + "0", + "-16128/390625", + "0", + "576/9765625" + ], + "('e', Fraction(1, 5))": [ + "-1196142592/375", + "0", + "-198682112/625", + "0", + "15364864/84375", + "0", + "-331776/78125", + "0", + "2304/390625" + ], + "('e', Fraction(3, 40))": [ + "1196142592/375", + "0", + "-1342950064/1875", + "0", + "28354646/84375", + "0", + "-4320861/312500", + "0", + "-5564463/12500000" + ], + "('eE', Fraction(1, 30), Fraction(1, 24))": [ + "0", + "-1059745792/1875", + "0", + "-4725213952/84375", + "0", + "402852608/6328125", + "0", + "-329472/390625", + "0", + "2304/1953125" + ], + "('eE', Fraction(1, 30), Fraction(1, 6))": [ + "0", + "-2119491584/1875", + "0", + "-9450427904/84375", + "0", + "805705216/6328125", + "0", + "-658944/390625", + "0", + "4608/1953125" + ], + "('eE', Fraction(1, 80), Fraction(1, 16))": [ + "0", + "794262672/625", + "0", + "-81791712/3125", + "0", + "24514137/156250", + "0", + "-681984/390625", + "0", + "-28117227/500000000" + ] + }, + "layer_cancellations": [ + [ + "1/24", + [ + "('1',)" + ] + ], + [ + "1/48", + [ + "('1',)", + "('E', Fraction(3, 16))", + "('K', Fraction(1, 48), Fraction(1, 24))", + "('K', Fraction(1, 48), Fraction(1, 6))" + ] + ], + [ + "1/80", + [ + "('1',)", + "('E', Fraction(1, 5))", + "('E', Fraction(3, 40))" + ] + ] + ], + "quadrature_remainder_sum": "6.21381774694680030565663976389602879124261835114676540685154E-55", + "quadrature_tail": "5.40488877623605478690138816214714996699398814519149302692619E-65" + }, + "owen_type_layer": { + "K_1/80_1/16": { + "1/4": { + "cumulative_layer": [ + "0.0312174775624259082976749478178245933483432339340057169299474", + "0.0312174775624259082976749478178245933483432339340057169299826" + ], + "owen_T": [ + "0.0312174775624259082976749478178245933483432339340057169121974", + "0.0312174775624259082976749478178245933483432339340057169474980" + ] + }, + "81/4": { + "cumulative_layer": [ + "20.5381371973009126057299000408000273773799756740867969013413", + "20.5381371973009126057299000408000273773799756740867969013649" + ], + "owen_T": [ + "20.5381371973009126057299000408000273773799756740867969012122", + "20.5381371973009126057299000408000273773799756740867969014984" + ] + }, + "100": { + "cumulative_layer": [ + "20.5765120396218307580906183508957531915225988992330229451570", + "20.5765120396218307580906183508957531915225988992330229452284" + ], + "owen_T": [ + "20.5765120396218307580906183508957531915225988992330229444124", + "20.5765120396218307580906183508957531915225988992330229461046" + ] + } + }, + "K_1/30_1/24": { + "1/4": { + "cumulative_layer": [ + "0.0312039330157050500612384331057801146378986915088572154703616", + "0.0312039330157050500612384331057801146378986915088572154703978" + ], + "owen_T": [ + "0.0312039330157050500612384331057801146378986915088572154561716", + "0.0312039330157050500612384331057801146378986915088572154842505" + ] + }, + "81/4": { + "cumulative_layer": [ + "11.2841167259181584699157774346793918427698323825872316848268", + "11.2841167259181584699157774346793918427698323825872316848415" + ], + "owen_T": [ + "11.2841167259181584699157774346793918427698323825872316847507", + "11.2841167259181584699157774346793918427698323825872316849183" + ] + }, + "100": { + "cumulative_layer": [ + "11.2841203268116122218023203328857390095151567735536801010267", + "11.2841203268116122218023203328857390095151567735536801010891" + ], + "owen_T": [ + "11.2841203268116122218023203328857390095151567735536801006328", + "11.2841203268116122218023203328857390095151567735536801014775" + ] + } + }, + "K_1/30_1/6": { + "1/4": { + "cumulative_layer": [ + "0.0311634049735495851207221365869627700973253813781908469259432", + "0.0311634049735495851207221365869627700973253813781908469259913" + ], + "owen_T": [ + "0.0311634049735495851207221365869627700973253813781908469190666", + "0.0311634049735495851207221365869627700973253813781908469328989" + ] + }, + "81/4": { + "cumulative_layer": [ + "7.71619021441145560725924210844156386474749302355847551756869", + "7.71619021441145560725924210844156386474749302355847551757486" + ], + "owen_T": [ + "7.71619021441145560725924210844156386474749302355847551751593", + "7.71619021441145560725924210844156386474749302355847551762407" + ] + }, + "100": { + "cumulative_layer": [ + "7.71619201485818653428398188158590744682097458721238360471594", + "7.71619201485818653428398188158590744682097458721238360472688" + ], + "owen_T": [ + "7.71619201485818653428398188158590744682097458721238360401738", + "7.71619201485818653428398188158590744682097458721238360544385" + ] + } + } + }, + "layers_interval_vs_float": { + "3": { + "F5_interval": [ + "5.7178385716165471595723471453580630280810165487E-8", + "5.7178385716165471595723471453580630280810404374E-8" + ], + "F5_direct_float_n20": "5.71783857161613e-08" + }, + "6": { + "F5_interval": [ + "0.4878428172838917427623082337202502226064969677430592229", + "0.4878428172838917427623082337202502226064969677434525232" + ], + "F5_direct_float_n20": "0.487842817283861" + }, + "10": { + "F5_interval": [ + "25055.3823953112841656522424296488692519216488248023079772594", + "25055.3823953112841656522424296488692519216488248023089773595" + ], + "F5_direct_float_n20": "25055.3823953095" + } + }, + "crosschecks": { + "m2 [0,1]": { + "T_quadrature": [ + "0.825283771446628026153734795989501098973937785745525336524934", + "0.825283771446628026153734795989501098973937785745525336565843" + ], + "Math197_closed_form": [ + "0.82528377144662802615373479598950109897393778574552533653260", + "0.82528377144662802615373479598950109897393778574552533653677" + ] + }, + "m2 [3,7/2]": { + "T_quadrature": [ + "0.0139921205104396079790627336205182811993011850118845207933045", + "0.0139921205104396079790627336205182811993011850118845208680820" + ], + "Math197_closed_form": [ + "0.0139921205104396079790627336205182811993011850118845208278966", + "0.0139921205104396079790627336205182811993011850118845208341167" + ] + }, + "m2 [-inf,0]": { + "T_quadrature": [ + "0.110138976347499685621780032029532612710583311728738766381446", + "0.110138976347499685621780032029532612710583311728738766385631" + ], + "Math197_closed_form": [ + "0.110138976347499685621780032029532612710583311728738766381762", + "0.110138976347499685621780032029532612710583311728738766382795" + ] + }, + "m3 full": [ + "5.32318026888989689492319873968726063578438852944092005011148", + "5.32318026888989689492319873968726063578438852944092005013162" + ], + "m3 F^(4)_[0,1]": [ + "0.153903300067681539201309170781277197985100298089699919851643", + "0.153903300067681539201309170781277197985100298089699919859783" + ], + "m4 full": [ + "14.2187634589773576013064656778740783173406023160820215841445", + "14.2187634589773576013064656778740783173406023160820215842808" + ], + "m4 F^(5)_[0,1]": [ + "0.0437174300836203305372424566264352813324007408534016628666462", + "0.0437174300836203305372424566264352813324007408534016628716563" + ] + }, + "float_control": { + "[0,1]": { + "float_eigenvalue_quadrature": "0.353131021", + "certified": [ + "0.353131021034538035163718545945115093462020092504131726835250", + "0.353131021034538035163718545945115093462020092504131727244403" + ] + }, + "[1,2]": { + "float_eigenvalue_quadrature": "8.980802568", + "certified": [ + "8.98080256786551439601323524155552010994814424212943227282558", + "8.98080256786551439601323524155552010994814424212943227406630" + ] + } + }, + "tilted": { + "v_mu_sigma2": [ + "17/25", + "33/34", + "133/5100" + ], + "eps": "0.108453837297795402480688369052653602949449849500832665110263", + "full_window": [ + "44.1306518750741327103620248954026331487822380630770420119935", + "44.1306518750741327103620248954026331487822380630770420128758" + ], + "windows": { + "[0,1]": { + "tilted_exact": [ + "0.458643108632784214791843630771388290132941986481399479425435", + "0.458643108632784214791843630771388290132941986481399479844088" + ], + "sandwich": [ + "0.0451885398213414111497882264588866799991655798173519896291277", + "2.68572653601198238492550903199083205666700009920450162509784" + ], + "reference": [ + "0.353131021034538035163718545945115093462020092504131726835250", + "0.353131021034538035163718545945115093462020092504131727244403" + ] + }, + "[1,2]": { + "tilted_exact": [ + "9.60771387239455380309599463055738068479526446373185971493242", + "9.60771387239455380309599463055738068479526446373185971617945" + ], + "sandwich": [ + "1.12778376443273512654947084830773764199132162400816340443417", + "68.3411120481063782323781448239473150456046653365943812593411" + ], + "reference": [ + "8.98080256786551439601323524155552010994814424212943227282558", + "8.98080256786551439601323524155552010994814424212943227406630" + ] + } + } + }, + "constants": { + "pi": [ + "3.14159265358979323846264338327950288419716939937510582097457", + "3.14159265358979323846264338327950288419716939937510582097528" + ], + "Gamma(7/6)": [ + "0.927719333630039200708349482534621018566466519145475576700351", + "0.927719333630039200708349482534621018566466519145475576939062" + ], + "(3/2)^(1/3)": [ + "1.14471424255333186780804221193967700891590692078793107209905", + "1.14471424255333186780804221193967700891590692078793107209906" + ] + }, + "exact_relations": { + "window coefficient": "c^(6)_{B,K} = 144 int_{S^5} p_G(0) p_V(0) I_B(u) J_K(u) dsigma(u), I_B = E[1{f in B} det(A)^2 1{A<0} | V=0] ([LP] (11.3), section 15)", + "reference": "(f, A) | V=0: f ~ N(0, 2/3), A = Q - f I, Q the 5x5 GOE; I^ref_B = (sqrt(3/8)/(46080 pi^(3/2))) int_0^inf e^{-3T^2/160} F_5(T) [Phi(2(b_+ - T/8)) - Phi(2(b_- - T/8))] dT", + "Owen-type kinds": "K_{al,g}(T) = int_0^T e^{-al w^2} E_g(w) dw = (pi/sqrt(al g)) [arctan(a)/(2 pi) - T_Owen(sqrt(2 al) T, a)], a = sqrt(g/al)", + "sandwich": "(1-eps) C'_ref <= C' <= (1+eps) C'_ref => (1-eps)^13/(1+eps)^8 I^ref_{B/sqrt(1-eps)} <= I'_B <= (1+eps)^13/(1-eps)^8 I^ref_{B/sqrt(1+eps)}", + "epsilon": "eps = ||C' - C'_ref||_F / lambda0, lambda0 = 23/100 <= lambda_min(C'_ref) = (6 - sqrt(92/3))/2", + "image bound": "E^(6)_L = 314496 (76 L^6 + 15) exp(-L^2/2)", + "c_6_ref": "Gamma(7/6) (3/2)^(1/3) D_5 / (64 sqrt3 pi^(7/2)) (Math-#204)" + } +} diff --git a/frontiers/c8_window_coefficient_d6_20261001/SOURCE_FILES.json b/frontiers/c8_window_coefficient_d6_20261001/SOURCE_FILES.json new file mode 100644 index 00000000..cc4d8cb4 --- /dev/null +++ b/frontiers/c8_window_coefficient_d6_20261001/SOURCE_FILES.json @@ -0,0 +1,137 @@ +{ + "schema": 1, + "object": "CL-C8-WINDOW-COEFFICIENT-D6-20261001-v1", + "scientific_effect": "NONE", + "certified": true, + "review_status": "AUTHOR_SIDE_CERTIFIED_ENCLOSURE_NONAUTHOR_READS_OPEN", + "disposition": "CERTIFIED_REFERENCE_AND_TORUS_ENCLOSURES_OF_THE_D6_WINDOW_COEFFICIENT_NO_REGISTER_CHANGE_C8_OPEN", + "declarative": true, + "executed": false, + "base": "124c37d9218120245098df453cf88d3e54b8c9ca", + "claim": "https://github.com/d6g8k5htny-coder/Math-/pull/213#issuecomment-5930520277", + "pins_on_main": [ + { + "label": "LP", + "path": "imports/lifetime_parent_20260925/UNIFORM_MATRIX_CAP_AND_LIFETIME.md", + "git_blob": "dfed3b8d318a3ab1950957f393307733a4bef3f2", + "role": "(11.3) the window coefficient in every d; section 15 factorization p_(f,V)(b,0) p_G(0) phi_tau(12k), z_0 = 36 k^2 E[det(A)^2 1{A<0} | f=b, V=0]", + "workflow_checked": true + }, + { + "label": "S24", + "path": "coefficients/side24_v1/PROOF.md", + "git_blob": "44b66f04f89fcd87383b3603fa69f1feb64cdddd", + "role": "sections 1-4: the reference law A = Q + sqrt(2/3) Z I, the image-bound pattern (76|x|^6 contraction bound, normalization), the covariance-sandwich scaling argument", + "workflow_checked": true + } + ], + "companion_not_on_main": [ + { + "pr": "https://github.com/d6g8k5htny-coder/Math-/pull/197", + "commit": "ebc93cd", + "path": "frontiers/c8_window_coefficient_planar_20260930/window_coefficient.py", + "git_blob": "43b4f51ba0e29c38468c40b4648f683848a3b218", + "role": "interval class, erf, arctan, Owen's T, incomplete gamma, gap factor G_K, library_exact (reused verbatim through Math-#205's copy; C41 PASS_TECHNICAL, review 5373681556); Owen's T also certifies the Owen-type kinds (rule OWEN_EXACT)" + }, + { + "pr": "https://github.com/d6g8k5htny-coder/Math-/pull/205", + "commit": "aebb4d4", + "path": "frontiers/c8_window_coefficient_d3_torus_20260930/window_coefficient_d3.py", + "git_blob": "ba24b3802eb0baf5ade65d49fd47ca862972f8bd", + "role": "the toolkit section and I3_ref (rule CROSSCHECK_M2_M4, m = 2) reused verbatim; Lemma S (its NOTE section 2) applied with n = 16, degree 10" + }, + { + "pr": "https://github.com/d6g8k5htny-coder/Math-/pull/199", + "commit": "617236d", + "path": "coefficients/side24_d4_20260930/cone_moment_d4.py", + "git_blob": "4787d6bf5a3bc469ae62075dde24e913b786e44a", + "role": "exact polynomial algebra (padd, pmul, pscale, pconst, pvar) reused verbatim" + }, + { + "pr": "https://github.com/d6g8k5htny-coder/Math-/pull/209", + "commit": "84a494a", + "path": "frontiers/c8_window_coefficient_d4_20261001/window_coefficient_d4.py", + "git_blob": "f77f7d2d158f9513f23f1b55fe0f839eb112d387", + "role": "the total-trace reduction and the Taylor quadrature with the Phi factor; the d = 4 band factor F^(4)_[0,1] (rule CROSSCHECK_M2_M4, m = 3)" + }, + { + "pr": "https://github.com/d6g8k5htny-coder/Math-/pull/213", + "commit": "b22ff4b", + "path": "frontiers/c8_window_coefficient_d5_20261001/window_coefficient_d5.py", + "git_blob": "da68e6ec22310c35904a70e203cc32cb9e5b19fc", + "role": "the m = 4 layers, generalized here to one layer recursion for every m <= 5; its stated F_4 (rule INNER_LAYERS) and the d = 5 band factor F^(5)_[0,1] (rule CROSSCHECK_M2_M4, m = 4); the tilted pair and the exact tests" + }, + { + "pr": "https://github.com/d6g8k5htny-coder/Math-/pull/204", + "commit": "5915af4", + "path": "coefficients/side24_d6_20260930/cone_moment_d6.py", + "git_blob": "c5e56694d96c885373a7b70868259eeb9fdfa746", + "role": "the certified cumulative layer for Owen-type kinds (method), Z_5 = 46080 pi^(3/2), the E_6 image-bound constant 314496" + }, + { + "pr": "https://github.com/d6g8k5htny-coder/Math-/pull/204", + "commit": "5915af4", + "path": "coefficients/side24_d6_20260930/RESULTS.json", + "git_blob": "031c8bbe38eeb473a2223cc58522803207b4d556", + "role": "the certified D_5 interval (rule D5_CERTIFIED) and the c_{6,ref} digits (rule C6_CONSISTENT)" + }, + { + "pr": "https://github.com/d6g8k5htny-coder/Math-/pull/202", + "commit": "4512309", + "path": "coefficients/side24_d5_closed_form_20260930/RESULTS.json", + "git_blob": "52ce9b5726b8ea92c429ea1193b1812869afe701", + "role": "the closed form of D_4 (rule CROSSCHECK_M2_M4, m = 4) and Z_4 = 1536 pi" + }, + { + "pr": "https://github.com/d6g8k5htny-coder/Math-/pull/201", + "commit": "3abcc10", + "path": "coefficients/side24_d4_closed_form_20260930/NOTE.md", + "git_blob": "cbf35dea", + "role": "the closed form of D_3 (rule CROSSCHECK_M2_M4, m = 3)" + } + ], + "files": [ + { + "path": "NOTE.md", + "bytes": 25762, + "sha256": "c143213807386240e612b6f3cc00b0f31dc6d2046f986b5f165dd2ab19422543", + "git_blob": "654c0aca5d04481418bd5768745100d338c9b544" + }, + { + "path": "RESULTS.json", + "bytes": 28348, + "sha256": "6d84046778b20102607550aa9763bd5e53ce38442f20721c33f86e5f9af07bbf", + "git_blob": "df361c43b691a90316b0f543150d1fd6935f897a" + }, + { + "path": "window_coefficient_d6.py", + "bytes": 70251, + "sha256": "a6c3717bc791f4f3b7fff171cf0c68426e757c684060a56ad4b59b5ea6737fe1", + "git_blob": "d38527e6490468c626072e299383f0652849a6bf" + } + ], + "revisions": [ + { + "version": "v1.1", + "finding": "band-share wording: xAI (Harper), this PR comment 5932429158; a scoped caption check on this PR, comment 5932380111 (its lane is not stated)", + "review": "https://github.com/d6g8k5htny-coder/Math-/pull/222#issuecomment-5932429158", + "change": "NOTE.md section 5: the column F^(d)_[0,1] relabelled as the height-window share (the same for every K); the C8 band's share of c_{d,ref}, F^(d)_[0,1] G_[1/2,2], added; section 10", + "unchanged": [ + "window_coefficient_d6.py", + "RESULTS.json", + ".github/workflows/c8-window-coefficient-d6.yml" + ] + }, + { + "version": "v1.2", + "finding": "W3 (Grok Bot agent 12): NOTE section 9 stale against section 10 and the PR body", + "review": "https://github.com/d6g8k5htny-coder/Math-/pull/222#issuecomment-5974639810", + "change": "NOTE.md section 9: completed review scopes and hosted runs point to the PR's disposition block; section 10", + "unchanged": [ + "window_coefficient_d6.py", + "RESULTS.json", + ".github/workflows/c8-window-coefficient-d6.yml" + ] + } + ] +} diff --git a/frontiers/c8_window_coefficient_d6_20261001/window_coefficient_d6.py b/frontiers/c8_window_coefficient_d6_20261001/window_coefficient_d6.py new file mode 100644 index 00000000..d38527e6 --- /dev/null +++ b/frontiers/c8_window_coefficient_d6_20261001/window_coefficient_d6.py @@ -0,0 +1,1591 @@ +"""The compact-window coefficient c_{B,K} of [LP] Theorem B in d = 6: reference value and transfer to the torus for every +L >= 10 and every frame. + +Standard library only. Run from anywhere: python -B -S window_coefficient_d6.py [--mutant NAME] (about two minutes) +Every printed interval [lo, hi] is a rigorous enclosure: interval arithmetic over `decimal` (60 digits, directed rounding; exp +and sqrt widened by two units in the last place and checked against exact rational brackets by rule LIBRARY_EXACT), pi by +Machin, fractional powers by exact integer roots, erf and the lower incomplete gamma function by positive series with geometric +tails (and, for |x| >= 12, the enclosure 0 < erfc(x) < e^{-x^2}/(x sqrt(pi))), arctan and Owen's T by alternating series bracketed +by consecutive partial sums (the Math-#197 toolkit, unchanged), and one one-dimensional integral by an order-48 Taylor expansion +at each cell centre with a Cauchy remainder (Math-#199 section 3), carrying the Owen-type kinds as a certified cumulative layer +(Math-#204 section 3). + +Object ([LP] = imports/lifetime_parent_20260925/UNIFORM_MATRIX_CAP_AND_LIFETIME.md, (11.3) with the section 15 factorization, +both stated for every d, for compact windows B = [b_-, b_+] and K = [k_-, k_+] of positive length, b_- < b_+ and +0 < k_- < k_+ < infinity; frame (u, w1, ..., w5), G = grad f, V = (f_uu, f_uw1, ..., f_uw5), A the 5 x 5 transverse Hessian, +t = f_uuu): + + c^(6)_{B,K} = 144 int_{S^5} p_G(0) p_V(0) I_B(u) J_K(u) dsigma(u), I_B(u) = E[1{f in B} det(A)^2 1{A < 0} | V = 0], + J_K(u) = int_K k^(4/3) phi_tau(12 k) dk, tau^2 = Var(t | G = 0). + +Reference kernel. (f, A) | V = 0 is f ~ N(0, 2/3) and A = Q - f I with Q the 5 x 5 GOE (density ~ exp(-tr Q^2/4)) independent +of f. On the ordered eigenvalue sector of A, with total trace -T and nested coordinates (w', w, s, q) (Jacobian 1/120), the +exponent is 3T^2/160 + w'^2/80 + w^2/48 + s^2/24 + q^2/8, the birth window enters only through +W_B(T) = Phi(2 (b_+ - T/8)) - Phi(2 (b_- - T/8)), and only odd powers of q occur, so the q-, s- and w-integrals are exact (in the +w-integral the Owen-type terms cancel identically, as in d = 5). The w'-integral leaves three Owen-type kinds +K_{al,g}(T) = int_0^T e^{-al w^2} E_g(w) dw, (al, g) = (1/80, 1/16), (1/30, 1/24), (1/30, 1/6), with E_g(w) = int_0^w e^{-g s^2} ds: + + I^ref_B = (sqrt(3/8)/(46080 pi^(3/2))) int_0^inf exp(-3T^2/160) F_5(T) W_B(T) dT, + F_5 = polynomials times exp(-T^2/5), exp(-3T^2/40), exp(-al T^2) E_g(T) and K_{al,g}(T), all three K with one polynomial R_10. + +The K_{al,g} are Owen's T functions, K = (pi/sqrt(al g)) [arctan(a)/(2 pi) - T(sqrt(2 al) T, a)], a = sqrt(g/al); the quadrature +carries them as a certified cumulative layer, and rule OWEN_EXACT checks that layer against Owen's series. B = R returns +Math-#204's D_5. The same code at m = 2 returns 29/6 - sqrt6 and Math-#197's d = 3 window values, at m = 3 Math-#201's D_3 and +Math-#209's d = 4 band factor, and at m = 4 Math-#202's D_4, Math-#213's stated F_4 and its d = 5 band factor. + +Torus. If (1-eps) C'_ref <= C' <= (1+eps) C'_ref for the covariance of the 16-vector (f, A) | V = 0, the density comparison and +the degree-10 homogeneity of det(A)^2 give (Lemma S of Math-#205 with n = 16, degree 10) + + (1-eps)^13/(1+eps)^8 I^ref_{B/sqrt(1-eps)} <= I'_B <= (1+eps)^13/(1-eps)^8 I^ref_{B/sqrt(1+eps)}. + +lambda_min(C'_ref) = (6 - sqrt(92/3))/2 = 0.2311 (below the d = 5 floor 1/4): eps = ||C' - C'_ref||_F / (23/100) over the entry +box of the d = 6 image bound E^(6)_L = 314496 (76 L^6 + 15) exp(-L^2/2); p_G(0), p_V(0) and tau^2 over the same box. One +enclosure covers every frame and every L >= 10. Windows touching k = 0 or infinity, unbounded B, and singleton windows +(b_- = b_+ or k_- = k_+, where the integral is zero) are evaluated as coefficient integrals only. +Scientific effect NONE. +""" +import argparse +import json +import math +import sys +import time +from decimal import Decimal, Context, ROUND_FLOOR, ROUND_CEILING, ROUND_HALF_EVEN, setcontext +from fractions import Fraction as Fr + +MUTANTS = ("image-shells", "d5-floor", "sandwich-swap", "window-scale", "trace-slope", "moment-recurrence", "remainder-dropped", + "parts-sign", "schur-sign", "tau-cross", "sphere", "tilt-gamma", "layer-dropped", "cumulative-shift") +MUT = None + +PREC = 60 +CF = Context(prec=PREC, rounding=ROUND_FLOOR, Emin=-9999999, Emax=9999999) +CC = Context(prec=PREC, rounding=ROUND_CEILING, Emin=-9999999, Emax=9999999) +CE = Context(prec=PREC, rounding=ROUND_HALF_EVEN, Emin=-9999999, Emax=9999999) +# Every arithmetic step below names one of the three contexts; operator syntax on Decimals would round to the thread's +# current context, which is widened as a safety net (rule LIBRARY_EXACT checks exp, sqrt and the negation). +setcontext(Context(prec=2 * PREC + 20, rounding=ROUND_HALF_EVEN, Emin=-9999999, Emax=9999999)) +ZERO, ONE, HALF = Decimal(0), Decimal(1), Decimal("0.5") +INF = math.inf +ROOT_DIGITS = 70 +K_ORDER, RHO, CELL, T_CUT = 48, Fr(6), Fr(1, 2), Fr(100) # Taylor order, Cauchy radius, cell width, cut (Math-#199 section 3) +D5_204 = ("44.1306518750741327103620248954026331485852349159356435082422", # Math-#204 RESULTS.json cone_moments/D_5 (certified) + "44.1306518750741327103620248954026331489792412102184405166029") +C6_REF_204 = "0.0076928786292368482786666312812273900811" # Math-#204 RESULTS.json (c_{6,ref}, certified there) +F4B_209 = "0.1539033000676815392013091707812771979851002980" # Math-#209 PINNED (F^(4)_[0,1], the d = 4 band birth factor) +F5B_213 = "0.04371743008362033053724245662643528133240074085" # Math-#213 PINNED (F^(5)_[0,1], the d = 5 band birth factor) +PINNED = {"c_6_ref": "0.007692878629236848278666631281227390081159238603", "D_5": "44.13065187507413271036202489540263314878223806307", "F6|[0,1]": "0.008001944363618010329693943918504838262453904498", "c6|[0,1]x[1/2,2] ref": "0.000004147268276047706822538464047034235965965259", "pi": "3.141592653589793238462643383279502884197169399375"} + + +# ---------------- interval arithmetic and special functions (Math-#197 toolkit, unchanged) +def ulp(x): + if x == 0: + return Decimal(1).scaleb(-PREC - 30) + return Decimal(1).scaleb(x.adjusted() - PREC + 1) + + +class Iv: + """Closed interval [lo, hi] of Decimals; every operation encloses the exact result.""" + __slots__ = ("lo", "hi") + + def __init__(self, lo, hi=None): + if hi is None: + hi = lo + if lo > hi: + raise ValueError("empty interval") + self.lo, self.hi = lo, hi + + @staticmethod + def frac(q): + q = Fr(q) + n, d = Decimal(q.numerator), Decimal(q.denominator) + return Iv(CF.divide(n, d), CC.divide(n, d)) + + @staticmethod + def of(x): + return x if isinstance(x, Iv) else Iv.frac(x) + + def __add__(a, b): + b = Iv.of(b) + return Iv(CF.add(a.lo, b.lo), CC.add(a.hi, b.hi)) + __radd__ = __add__ + + def __neg__(a): + return Iv(a.hi.copy_negate(), a.lo.copy_negate()) # exact: unary minus would round to the default context + + def __sub__(a, b): + b = Iv.of(b) + return Iv(CF.subtract(a.lo, b.hi), CC.subtract(a.hi, b.lo)) + + def __rsub__(a, b): + return Iv.of(b) - a + + def __mul__(a, b): + b = Iv.of(b) + c = [(x, y) for x in (a.lo, a.hi) for y in (b.lo, b.hi)] + return Iv(min(CF.multiply(x, y) for x, y in c), max(CC.multiply(x, y) for x, y in c)) + __rmul__ = __mul__ + + def __truediv__(a, b): + b = Iv.of(b) + if b.lo <= 0 <= b.hi: + raise ZeroDivisionError("interval division by an interval containing 0") + c = [(x, y) for x in (a.lo, a.hi) for y in (b.lo, b.hi)] + return Iv(min(CF.divide(x, y) for x, y in c), max(CC.divide(x, y) for x, y in c)) + + def __rtruediv__(a, b): + return Iv.of(b) / a + + def sq(a): + if a.lo <= 0 <= a.hi: + return Iv(ZERO, max(CC.multiply(a.lo, a.lo), CC.multiply(a.hi, a.hi))) + return a * a + + def exp(a): + lo, hi = CE.exp(a.lo), CE.exp(a.hi) + return Iv(max(ZERO, CF.subtract(lo, CC.multiply(Decimal(2), ulp(lo)))), CC.add(hi, CC.multiply(Decimal(2), ulp(hi)))) + + def sqrt(a): + if a.lo < 0: + raise ValueError("sqrt of an interval with negative part") + lo, hi = CE.sqrt(a.lo), CE.sqrt(a.hi) + return Iv(max(ZERO, CF.subtract(lo, CC.multiply(Decimal(2), ulp(lo)))), CC.add(hi, CC.multiply(Decimal(2), ulp(hi)))) + + def width(a): + return CC.subtract(a.hi, a.lo) + + def contains(a, x): + x = Iv.of(x) + return a.lo <= x.lo and x.hi <= a.hi + + def intersects(a, b): + return a.lo <= b.hi and b.lo <= a.hi + + def mag(a): + return max(abs(a.lo), abs(a.hi)) + + def mid_float(a): + return float((a.lo + a.hi) / 2) + + def pair(a): + return [str(a.lo), str(a.hi)] + + +def hull(a, b): + return Iv(min(a.lo, b.lo), max(a.hi, b.hi)) + + +def inroot(N, n): + """floor(N ** (1/n)) for integers N >= 0, n >= 1: integer Newton iteration started above the root.""" + if N < 2: + return N + x = 1 << ((N.bit_length() + n - 1) // n) + while True: + y = ((n - 1) * x + N // pow(x, n - 1)) // n + if y >= x: + break + x = y + if not (pow(x, n) <= N < pow(x + 1, n)): + raise ValueError("integer root failed") + return x + + +def root_frac(q, n): + """Rigorous bracket of q ** (1/n) for a positive rational q, width 10^-ROOT_DIGITS.""" + q = Fr(q) + if q <= 0: + raise ValueError("root of a nonpositive number") + p, d = q.numerator, q.denominator + scale = 10 ** ROOT_DIGITS + r = inroot(p * d ** (n - 1) * scale ** n, n) + den = d * scale + return Iv(Iv.frac(Fr(r, den)).lo, Iv.frac(Fr(r + 1, den)).hi) + + +def root_iv(x, n): + """x ** (1/n) for an interval x > 0 (monotone: endpoint brackets).""" + return Iv(root_frac(Fr(x.lo), n).lo, root_frac(Fr(x.hi), n).hi) + + +def pi_interval(): + def atan_inv(n, terms): + x = Fr(1, n) + partial, sums = Iv(ZERO), [] + for m in range(terms + 2): + t = Iv.frac(x ** (2 * m + 1) / (2 * m + 1)) + partial = partial + t if m % 2 == 0 else partial - t + sums.append(partial) + return hull(sums[-1], sums[-2]) # alternating series with decreasing terms + return atan_inv(5, 44) * 16 - atan_inv(239, 22) * 4 + + +PI = pi_interval() +SQRT2 = Iv.frac(2).sqrt() +SQRT3 = Iv.frac(3).sqrt() +SQRT2PI = (2 * PI).sqrt() +SQRTPI = PI.sqrt() +TINY = Decimal(10) ** (-PREC - 8) + + +def erf_point(x, coarse=False): + """erf(x) at a Decimal point: (2/sqrt pi) exp(-x^2) sum_n 2^n x^(2n+1)/(2n+1)!!, positive terms; once the term ratio + 2x^2/(2n+3) is at most 1/2 the omitted tail is at most twice the first omitted term. coarse=True stops at that point + (a valid, wider enclosure used by TRUNCATION_NESTING); the production run continues until the tail is negligible.""" + tail = MUT != "tail-dropped" + coarse = coarse or not tail + if x < 0: + return -erf_point(-x, coarse) + if x == 0: + return Iv(ZERO) + X = Iv(x) + x2 = X * X + t, S, n = X, Iv(ZERO), 0 + while True: + S = S + t + ratio = x2 * 2 / (2 * n + 3) + t = t * ratio + n += 1 + if ratio.hi <= HALF and (coarse or t.hi < S.lo * TINY): + break + if n > 200000: + raise RuntimeError("erf series did not converge") + rest = Iv(ZERO, CC.multiply(Decimal(2), t.hi)) if tail else Iv(ZERO) + return (S + rest) * (-x2).exp() * 2 / SQRTPI + + +def erf_iv(x, coarse=False): + return Iv(erf_point(x.lo, coarse).lo, erf_point(x.hi, coarse).hi) + + +def Phi(x): + return (1 + erf_iv(x / SQRT2)) / 2 + + +def phi(x): + return (-(x.sq()) / 2).exp() / SQRT2PI + + +def atan_iv(a, N=90): + """arctan(a) for an interval 0 <= a <= 1: halve the argument, then the alternating series bracketed by consecutive + partial sums (for each fixed argument the terms decrease, so the value lies between them).""" + a1 = a / (1 + (1 + a.sq()).sqrt()) + a1sq, p, S, prev = a1.sq(), a1, Iv(ZERO), None + for n in range(N + 1): + t = p / (2 * n + 1) + prev = S + S = S + t if n % 2 == 0 else S - t + p = p * a1sq + return 2 * hull(prev, S) + + +def owen_T(h, a, N=160, bracket=True): + """Owen's T(h, a) = (1/2 pi) int_0^a exp(-h^2 (1+x^2)/2)/(1+x^2) dx for 0 <= a <= 1, by Owen's series + T = (1/2 pi)[arctan a - sum_j c_j a^(2j+1)], c_j = (-1)^j/(2j+1) [1 - exp(-h^2/2) sum_{i<=j} (h^2/2)^i/i!]. + The bracket is P(Poisson(h^2/2) > j), decreasing in j, so for fixed (h, a) the terms alternate with decreasing size and + the sum lies between consecutive partial sums; interval evaluation of both encloses it for every point of the intervals.""" + if MUT == "no-owen": + return Iv(ZERO) + h2 = h.sq() / 2 + e = (-h2).exp() + pois, term, apow, a2 = Iv(ONE), Iv(ONE), a, a.sq() + S, prev = Iv(ZERO), None + for j in range(N + 1): + if j > 0: + term = term * h2 / j + pois = pois + term + cj = apow / (2 * j + 1) * (1 - e * pois) + prev = S + S = S + cj if j % 2 == 0 else S - cj + apow = apow * a2 + val = hull(prev, S) if bracket else S + return (atan_iv(a) - val) / (2 * PI) + + +def owen_T_eval(h, a): + if MUT == "tail-dropped": + return owen_T(h, a, N=6, bracket=False) + return owen_T(h, a) + + +def SHAPE(): + """Shape parameter of the gap-window incomplete gamma factor (mutant gamma-shape perturbs it).""" + return Fr(4, 3) if MUT == "gamma-shape" else Fr(7, 6) + + +def gamma_lower_point(x, a, coarse=False): + """gamma(a, x) = int_0^x s^(a-1) e^-s ds = x^a e^-x sum_n x^n/(a (a+1) ... (a+n)), positive terms; once the term ratio + x/(a+n) is at most 1/2 the omitted tail is at most twice the first omitted term. x a Decimal point >= 0, a = 7/6 or 4/3.""" + tail = MUT != "tail-dropped" + coarse = coarse or not tail + if x == 0: + return Iv(ZERO) + X = Iv(x) + t, S, n = Iv.frac(1 / a), Iv(ZERO), 0 + while True: + S = S + t + n += 1 + ratio = X / Iv.frac(a + n) + t = t * ratio + if ratio.hi <= HALF and (coarse or t.hi < S.lo * TINY): + break + if n > 200000: + raise RuntimeError("incomplete gamma series did not converge") + rest = Iv(ZERO, CC.multiply(Decimal(2), t.hi)) if tail else Iv(ZERO) + xa = X * root_frac(Fr(x), 6) if a == Fr(7, 6) else X * root_frac(Fr(x), 3) # x^a + return xa * (-X).exp() * (S + rest) + + +def gamma_lower(x, a, coarse=False): + """gamma(a, x) for an interval x >= 0 or x = INF (monotone in x).""" + if isinstance(x, float): + return GAMMA[a] + return Iv(gamma_lower_point(x.lo, a, coarse).lo, gamma_lower_point(x.hi, a, coarse).hi) + + +def gamma_complete(a, X=Decimal(120)): + """Gamma(a) = gamma(a, X) + Gamma(a, X), 1 < a < 2, with X^(a-1) e^-X <= Gamma(a, X) <= X^(a-1) e^-X / (1 - (a-1)/X).""" + low = gamma_lower_point(X, a) + xa1 = root_frac(Fr(X), 6) if a == Fr(7, 6) else root_frac(Fr(X), 3) + up = xa1 * (-Iv(X)).exp() + return low + Iv(up.lo, (up / (1 - Iv.frac(a - 1) / Iv(X))).hi) + + +GAMMA = {Fr(7, 6): gamma_complete(Fr(7, 6)), Fr(4, 3): gamma_complete(Fr(4, 3))} +GAMMA76 = GAMMA[Fr(7, 6)] +CBRT32 = root_frac(Fr(12), 3) / 2 # (3/2)^(1/3) = 12^(1/3)/2 +C72 = 1 / (72 * root_frac(Fr(72), 6)) # 72^(-7/6) +SIGMA_REF = Iv.frac(Fr(2, 3)).sqrt() # sd of f | V=0 for the reference kernel +D2_EXACT = Iv.frac(Fr(29, 6)) - Iv.frac(6).sqrt() # SIDE24 cone moment m = 2 +C3_REF = GAMMA76 * CBRT32 * D2_EXACT / (2 * SQRT3 * PI * PI * SQRTPI) # SIDE24 (1), d = 3 + + +# ---------------- the Gaussian integrals int_{-inf}^x t^n phi(t) Phi(a t) dt, n <= 4 (Math-#197 section 2, unchanged) +def G_finite(x, a): + rho2 = 1 + a.sq() + rho = rho2.sqrt() + Phx, phx, Phax = Phi(x), phi(x), Phi(a * x) + T = owen_T_eval(x, a) + Erho = (-(rho2 * x.sq()) / 2).exp() + Phrx = Phi(rho * x) + x2 = x.sq() + G0 = Phx / 2 - T + G1 = -phx * Phax + a / (rho * SQRT2PI) * Phrx + G2 = G0 - x * phx * Phax - a / (2 * PI * rho2) * Erho + G3 = -x2 * phx * Phax + a / (2 * PI) * (-(x / rho2) * Erho + SQRT2PI / (rho2 * rho) * Phrx) + 2 * G1 + G4 = -x * x2 * phx * Phax - a / (2 * PI) * (x2 / rho2 + 2 / rho2.sq()) * Erho + 3 * G2 + return [G0, G1, G2, G3, G4] + + +def G_limit(a, sign): + if sign < 0: + return [Iv(ZERO)] * 5 + rho2 = 1 + a.sq() + rho = rho2.sqrt() + g1 = a / (rho * SQRT2PI) + return [Iv.frac(Fr(1, 2)), g1, Iv.frac(Fr(1, 2)), a / (rho2 * rho * SQRT2PI) + 2 * g1, Iv.frac(Fr(3, 2))] + + +def G_at(x, a): + """x an interval or +-INF.""" + if isinstance(x, float): + return G_limit(a, 1 if x > 0 else -1) + return G_finite(x, a) + + +def P_reg(s): + """Regularized lower incomplete gamma P(7/6, s) (s an interval >= 0 or INF).""" + a = SHAPE() + return gamma_lower(s, a) / GAMMA76 + + +def G_ref(K): + lo, hi = K + def s(k): + return k if isinstance(k, float) else 12 * Iv.frac(k * k) + return P_reg(s(hi)) - P_reg(s(lo)) + + +def library_exact(): + """exp and sqrt of the decimal module, and the interval negation, against exact rational brackets: e^-q by the + rational series of e^(-q/64) raised to the 64th power (alternating tail below the last term), sqrt by an integer root.""" + ok = True + args = [Fr(48) - Fr(51035039649098175696, 10 ** 29), Fr(3), Fr(120), Fr(1, 2), Fr(288), Fr(5, 7)] + for q in args: + g = -q / 64 + parts = [] + term, ssum = Fr(1), Fr(0) + for k in range(0, 60): + ssum += term + parts.append(ssum) + term = term * g / (k + 1) + lo, hi = min(parts[-1], parts[-2]), max(parts[-1], parts[-2]) # alternating series with decreasing terms + exact = Iv(Iv.frac(lo ** 64).lo, Iv.frac(hi ** 64).hi) # e^-q = (e^(-q/64))^64, both endpoints positive + X = Iv.frac(q) + got = (-X).exp() + ok &= got.intersects(exact) and CC.divide(got.width(), got.lo) < Decimal(10) ** (-PREC + 3) + r = root_frac(q, 2) + ok &= X.sqrt().intersects(r) and CC.divide(X.sqrt().width(), r.lo) < Decimal(10) ** (-PREC + 3) + ok &= (-X).lo == X.hi.copy_negate() and (-X).hi == X.lo.copy_negate() + return ok + + +# ---------------- the d = 3 reference birth integral with interval endpoints (Math-#197 (4.1) and its three primitives) +def endpoint(b): + return b if isinstance(b, (float, Iv)) else Iv.frac(b) + + +def I3_ref(lo, hi): + """int_lo^hi phi_{2/3}(b) m_{3,b} db, m_{3,b} = (b^3 + b) phi(b) + (b^4 + 2 b^2 + 7) Phi(b) - 4 sqrt2 exp(-b^2/4) Phi(b/sqrt2); + lo, hi are rationals, intervals or +-INF.""" + sig = SIGMA_REF + def term_i(b): + if isinstance(b, float): + return Iv(ZERO) + return -(b.sq() * 2 / 5 + Iv.frac(Fr(18, 25))) * (-(b.sq()) * 5 / 4).exp() / (2 * PI * sig) + def term_ii(b): + G = G_at(b if isinstance(b, float) else b / sig, sig) + return sig.sq().sq() * G[4] + 2 * sig.sq() * G[2] + 7 * G[0] + def term_iii(b): + G = G_at(b if isinstance(b, float) else b * SQRT2, Iv.frac(Fr(1, 2))) + return 4 / sig * G[0] + lo, hi = endpoint(lo), endpoint(hi) + cross = Iv.frac(Fr(1, 2)) if MUT == "d3-cross" else Iv(ONE) + return (term_i(hi) - term_i(lo)) + (term_ii(hi) - term_ii(lo)) - cross * (term_iii(hi) - term_iii(lo)) + + +# ---------------- exact polynomial algebra (Math-#199 cone_moment_d4.py, unchanged) +def padd(u, v): + w = dict(u) + for k, c in v.items(): + w[k] = w.get(k, Fr(0)) + c + return {k: c for k, c in w.items() if c != 0} + + +def pmul(u, v): + w = {} + for k1, c1 in u.items(): + for k2, c2 in v.items(): + k = tuple(a + b for a, b in zip(k1, k2)) + w[k] = w.get(k, Fr(0)) + c1 * c2 + return {k: c for k, c in w.items() if c != 0} + + +def pscale(u, c): + return {k: v * c for k, v in u.items() if v * c != 0} + + +def pconst(c, nvars): + return {(0,) * nvars: Fr(c)} if c != 0 else {} + + +def pvar(i, nvars): + e = [0] * nvars + e[i] = 1 + return {tuple(e): Fr(1)} + + +def plin(coeffs): + """sum_i coeffs[i] * variable_i.""" + out = {} + for i, c in enumerate(coeffs): + if c: + out = padd(out, pscale(pvar(i, len(coeffs)), Fr(c))) + return out + + +def lpadd(u, v): + n = max(len(u), len(v)) + return [(u[i] if i < len(u) else Fr(0)) + (v[i] if i < len(v) else Fr(0)) for i in range(n)] + + +def lpmul(u, v): + if not u or not v: + return [] + w = [Fr(0)] * (len(u) + len(v) - 1) + for i, a in enumerate(u): + for j, b in enumerate(v): + w[i + j] += a * b + return w + + +# ---------------- the reference birth integral as a one-dimensional integral in the total trace (NOTE section 2) +def Jk_closed(k): + """J_k(s) = int_0^s q^k e^{-q^2/8} dq for odd k: J_k = B_k + C_k(s) e^{-s^2/8} with B_k rational and C_k a rational polynomial + (list by power): J_1 = 4 - 4 e^{-s^2/8}, J_{k+2} = -4 s^{k+1} e^{-s^2/8} + 4 (k+1) J_k.""" + if k % 2 != 1: + raise ValueError("even power of q") + B, C, j = Fr(4), [Fr(-4)], 1 + while j < k: + C2 = [4 * (j + 1) * c for c in C] + [Fr(0)] * (j + 2 - len(C)) + C2[j + 1] += -4 + B, C, j = 4 * (j + 1) * B, C2, j + 2 + return B, C + + +def moments(gamma, nmax): + """M_j(T) = int_0^T s^j e^{-gamma s^2} ds = alpha_j + rho_j(T) e^{-gamma T^2} + beta_j E_gamma(T), E_gamma = (1/2) sqrt(pi/gamma) + erf(sqrt(gamma) T): M_0 = E_gamma, M_1 = (1 - e^{-gamma T^2})/(2 gamma), M_{j+2} = ((j+1) M_j - T^{j+1} e^{-gamma T^2})/(2 gamma).""" + g = Fr(gamma) + M = [(Fr(0), [], Fr(1)), (1 / (2 * g), [-1 / (2 * g)], Fr(0))] + for j in range(nmax - 1): + a, r, b = M[j] + f = Fr(j if MUT == "moment-recurrence" else j + 1) + r2 = [x * f / (2 * g) for x in r] + [Fr(0)] * (j + 2 - len(r)) + r2[j + 1] += -1 / (2 * g) + M.append((f * a / (2 * g), r2, f * b / (2 * g))) + return M + + +F_MEMO = {} + + +def F_pieces(m): + if m not in F_MEMO: + F_MEMO[m] = F_pieces_build(m) + return F_MEMO[m] + + +def F_pieces_build(m): + """The inner integrals in closed form, for every m = 2, ..., 5 by one generic layer recursion (NOTE section 2). On the ordered + eigenvalue sector mu_m = -a, mu_{m-1} = -a - p_1, ..., mu_1 = -a - p_1 - ... - p_{m-1} put x_0 = T = m a + (m-1) p_1 + ... + p_{m-1} + (minus the trace) and x_k = (m-k) p_k + ... + p_{m-1}, k = 1, ..., m-1 (Jacobian 1/m!, region T > x_1 > ... > x_{m-1} > 0); the + integrand is P(x) e^{-sum_k x_k^2 c_k} with c_{m-1} = 1/8 and c_{m-2}, c_{m-3}, c_{m-4} = 1/24, 1/48, 1/80 (eight times the + triangular numbers), P = prod mu_i^2 times the Vandermonde. The q = x_{m-1} layer is exact + (Math-#201 L1); every further layer integrates y^j e^{-c y^2} times a kind ('1', ('e', g), ('E', g), ('eE', al, g)) from 0 to + the next variable by moment_dict and L_dict. An Owen-type kind ('K', al, g) entering a further layer is refused (ValueError). + Returns ({kind: polynomial in T, list by power}, outer gamma 3/(4 m (m+3)), [(c, generated kinds that came out zero) per layer]).""" + nv = m + xs = [pvar(i, nv) for i in range(nv)] + + def lin(i, j, den): + co = [Fr(0)] * nv + co[i] += Fr(1, den) + if j is not None: + co[j] -= Fr(1, den) + return plin(co) + a = lin(0, 1, m) + ps = [lin(k, k + 1 if k + 1 < m else None, m - k) for k in range(1, m)] + mus, run = [pscale(a, Fr(-1))], a + for p in ps: + run = padd(run, p) + mus.append(pscale(run, Fr(-1))) + P = {(0,) * nv: Fr(1)} + for mu in mus: + P = pmul(P, pmul(mu, mu)) + for i in range(m): + for j in range(i + 1, m): + P = pmul(P, padd(mus[i], pscale(mus[j], Fr(-1)))) # mus[i] > mus[j] + U, V = {}, {} + for mono, c in P.items(): + Bk, Ck = Jk_closed(mono[-1]) # raises on an even power of q + base = mono[:-1] + U[base] = U.get(base, Fr(0)) + c * Bk + for n, x in enumerate(Ck): + if x: + mm = list(base) + mm[-1] += n + mm = tuple(mm) + V[mm] = V.get(mm, Fr(0)) + c * x + pieces = {("1",): {k: v for k, v in U.items() if v}, ("e", Fr(1, 8)): {k: v for k, v in V.items() if v}} + stages, nvc = [], nv - 1 + for c in (Fr(1, 24), Fr(1, 48), Fr(1, 80))[:m - 2]: + pieces = layer(pieces, nvc, c) + nvc -= 1 + if MUT == "layer-dropped" and c == Fr(1, 48): + pieces.pop(("E", Fr(3, 16)), None) # a cancelled layer silently not derived + stages.append((c, sorted(("%s" % (k,)) for k, p in pieces.items() if not p))) + out = {} + for kind, poly in pieces.items(): + lst = [] + for mono, v in poly.items(): + n = mono[0] + lst += [Fr(0)] * (n + 1 - len(lst)) + lst[n] += v + out[kind] = lst + return out, Fr(3, 4 * m * (m + 3)), stages + + +def layer(pieces, nv, c): + """Integrate the last variable y of nv variables (T, x_1, ..., y) over (0, x_{nv-2}) against e^{-c y^2}, kind by kind.""" + out = {} + for kind, poly in pieces.items(): + for mono, coef in poly.items(): + j = mono[-1] + if kind[0] == "1": + lay = moment_dict(c, j) + elif kind[0] == "e": + lay = moment_dict(c + kind[1], j) + elif kind[0] == "E": + lay = L_dict(c, kind[1], j) + elif kind[0] == "eE": + lay = L_dict(c + kind[1], kind[2], j) + else: + raise ValueError("an Owen-type kind %r survived an inner layer" % (kind,)) + base = mono[:-1] + for k2, pl in lay.items(): + acc = out.setdefault(k2, {}) + for n, x in enumerate(pl): + if x: + mm = list(base) + mm[-1] += n + mm = tuple(mm) + acc[mm] = acc.get(mm, Fr(0)) + coef * x + return {k: {mm: v for mm, v in p.items() if v} for k, p in out.items()} + + +def moment_dict(gamma, j): + """M_j(gamma, T) = int_0^T w^j e^{-gamma w^2} dw by kind: ('1',) constant, ('e', gamma) e^{-gamma T^2}, ('E', gamma) E_gamma(T).""" + a_, r_, b_ = moments(gamma, max(j, 1))[j] + out = {} + if a_: + out[("1",)] = [a_] + if any(r_): + out[("e", gamma)] = list(r_) + if b_: + out[("E", gamma)] = [b_] + return out + + +L_MEMO = {} + + +def L_dict(al, g, j): + """L_j(T) = int_0^T w^j e^{-al w^2} E_g(w) dw, E_g(w) = (1/2) sqrt(pi/g) erf(sqrt(g) w): L_0 = Kfun_{al,g}(T) and, by parts, + L_j = -(1/(2 al)) T^{j-1} e^{-al T^2} E_g(T) + ((j-1)/(2 al)) L_{j-2} + (1/(2 al)) M_{j-1}(al + g, T) + (the last coefficient is (1/2) sqrt(pi/g) . sqrt(g)/(al sqrt(pi)) = 1/(2 al), rational).""" + key = (al, g, j) + if key in L_MEMO: + return L_MEMO[key] + if j == 0: + res = {("K", al, g): [Fr(1)]} + else: + res = {} + + def add(kind, poly): + res[kind] = lpadd(res.get(kind, []), poly) + add(("eE", al, g), [Fr(0)] * (j - 1) + [(1 if MUT == "parts-sign" else -1) / (2 * al)]) + if j >= 2: + for kind, pl in L_dict(al, g, j - 2).items(): + add(kind, [x * Fr(j - 1) / (2 * al) for x in pl]) + for kind, pl in moment_dict(al + g, j - 1).items(): + add(kind, [x / (2 * al) for x in pl]) + L_MEMO[key] = res + return res + + +def trim(poly): + poly = list(poly) + while poly and poly[-1] == 0: + poly.pop() + return poly + + +def g_terms(m, gamma_out=None): + """g_m(T) = e^{-gamma_out T^2} F_m(T) as terms (const, polynomial, gamma, kappa, K): const * poly(T) e^{-gamma T^2} [erf(kappa T)] + [K_{al,g}(T)], with K_{al,g}(T) = int_0^T e^{-al w^2} E_g(w) dw the Owen-type kinds (certified cumulative layer, NOTE section 3).""" + pieces, g_ref, _ = F_pieces(m) + gout = g_ref if gamma_out is None else gamma_out + terms = [] + for key, poly in sorted(pieces.items(), key=lambda kv: str(kv[0])): + poly = trim(poly) + if not poly: + continue + if key[0] == "1": + terms.append((Iv(ONE), poly, Iv.frac(gout), None, None)) + elif key[0] == "e": + terms.append((Iv(ONE), poly, Iv.frac(gout + key[1]), None, None)) + elif key[0] == "E": + terms.append(((PI / Iv.frac(key[1])).sqrt() / 2, poly, Iv.frac(gout), Iv.frac(key[1]).sqrt(), None)) + elif key[0] == "eE": + terms.append(((PI / Iv.frac(key[2])).sqrt() / 2, poly, Iv.frac(gout + key[1]), Iv.frac(key[2]).sqrt(), None)) + else: + terms.append((Iv(ONE), poly, Iv.frac(gout), None, (key[1], key[2]))) + return terms + + +def K_float(al, g, T, n=4000): + """Float control of K_{al,g}(T) by Simpson's rule.""" + if T == 0: + return 0.0 + h = T / n + f = [math.exp(-al * (i * h) ** 2) * 0.5 * math.sqrt(math.pi / g) * math.erf(math.sqrt(g) * i * h) for i in range(n + 1)] + return h / 3 * (f[0] + f[-1] + 4 * sum(f[1:-1:2]) + 2 * sum(f[2:-1:2])) + + +def F_float(pieces, T): + tot = 0.0 + for key, poly in pieces.items(): + pv = sum(float(c) * T ** i for i, c in enumerate(poly)) + if key[0] == "1": + tot += pv + elif key[0] == "e": + tot += pv * math.exp(-float(key[1]) * T * T) + elif key[0] == "E": + g = float(key[1]) + tot += pv * 0.5 * math.sqrt(math.pi / g) * math.erf(math.sqrt(g) * T) + elif key[0] == "eE": + al, g = float(key[1]), float(key[2]) + tot += pv * math.exp(-al * T * T) * 0.5 * math.sqrt(math.pi / g) * math.erf(math.sqrt(g) * T) + elif any(poly): + tot += pv * K_float(float(key[1]), float(key[2]), T) + return tot + + +# ---------------- Taylor quadrature with Cauchy remainders (Math-#199 section 3; factors as functions, plus Phi(alpha - beta T)) +def ipow(x, n): + r = Iv(ONE) + for _ in range(n): + r = r * x + return r + + +ERF_FAST = Decimal(12) + + +def erf_point_fast(x): + """erf at a Decimal point: the toolkit series for |x| < 12; for |x| >= 12 the enclosure [1 - e^{-x^2}/(x sqrt(pi)), 1] (with sign), + since 0 < erfc(x) < e^{-x^2}/(x sqrt(pi)) for x > 0 (here below 1.5e-64).""" + if abs(x) < ERF_FAST: + return erf_point(x) + X = Iv(abs(x)) + t = (-(X.sq())).exp() / (X * SQRTPI) + v = Iv((1 - t).lo, ONE) + return v if x > 0 else -v + + +def erf_fast(x): + return Iv(erf_point_fast(x.lo).lo, erf_point_fast(x.hi).hi) + + +def Phi_fast(x): + return (1 + erf_fast(x / SQRT2)) / 2 + + +CACHE = {} + + +def gauss_coeffs(g, s0, K): + """Taylor coefficients of exp(-g s^2) at s0: (n+1) y_{n+1} = -2 g (s0 y_n + y_{n-1}).""" + key = ("g", g.lo, g.hi, s0, K) + if key in CACHE: + return CACHE[key] + S0 = Iv.frac(s0) + y = [(-(g * S0.sq())).exp()] + prev = Iv(ZERO) + for n in range(K): + nxt = -(2 * g) * (S0 * y[n] + prev) / (n + 1) + prev = y[n] + y.append(nxt) + CACHE[key] = y + return y + + +def gauss_bound(g, s0, rho): + inner = max(Fr(0), abs(Fr(s0)) - Fr(rho)) + return (g * Iv.frac(rho * rho) - g * Iv.frac(inner * inner)).exp() # |e^{-g z^2}| = e^{-g(x^2 - y^2)}, |y| <= rho + + +def erf_coeffs(kappa, s0, K): + """erf(kappa s): value at s0, then (2 kappa/sqrt pi) times the Gaussian's coefficients divided by n + 1.""" + key = ("e", kappa.lo, kappa.hi, s0, K) + if key in CACHE: + return CACHE[key] + y = gauss_coeffs(kappa.sq(), s0, K) + pref = 2 * kappa / SQRTPI + out = [erf_fast(kappa * Iv.frac(s0))] + [pref * y[n] / (n + 1) for n in range(K)] + CACHE[key] = out + return out + + +def erf_bound(kappa, s0, rho): + return 2 * kappa / SQRTPI * Iv.frac(abs(Fr(s0)) + Fr(rho)) * (kappa.sq() * Iv.frac(Fr(rho) ** 2)).exp() + + +def poly_coeffs(c, s0, K): + out = [] + for n in range(K + 1): + v = Fr(0) + for j, cj in enumerate(c): + if j >= n: + v += cj * math.comb(j, n) * Fr(s0) ** (j - n) + out.append(Iv.frac(v)) + return out + + +def poly_bound(c, s0, rho): + r = abs(Fr(s0)) + Fr(rho) + return Iv.frac(sum(abs(cj) * r ** j for j, cj in enumerate(c))) + + +def t_mul(u, v): + K = len(u) - 1 + out = [] + for n in range(K + 1): + acc = Iv(ZERO) + for i in range(n + 1): + acc = acc + u[i] * v[n - i] + out.append(acc) + return out + + +def phi_aff_coeffs(alpha, beta, s0, K): + """Phi(alpha - beta T) at T = s0 + x (alpha, beta real intervals, beta > 0): c_0 = Phi(w0), w0 = alpha - beta s0, and + c_{n+1} = -(beta/sqrt(2 pi)) y_n/(n+1) with y the coefficients of exp(-(w0 - beta x)^2/2), (n+1) y_{n+1} = beta w0 y_n - beta^2 y_{n-1}.""" + w0 = alpha - beta * Iv.frac(s0) + y = [(-(w0.sq()) / 2).exp()] + prev = Iv(ZERO) + b2 = beta.sq() + for n in range(K): + nxt = (beta * w0 * y[n] - b2 * prev) / (n + 1) + prev = y[n] + y.append(nxt) + pref = beta / SQRT2PI + return [Phi_fast(w0)] + [-(pref * y[n] / (n + 1)) for n in range(K)] + + +def phi_aff_bound(beta, rho): + """|Phi(alpha - beta z)| <= Phi(w0) + |int_{w0}^{w} phi| <= 1 + beta rho e^{beta^2 rho^2/2}/sqrt(2 pi) on |z - s0| <= rho + (the segment from the real point w0 has |Im t| <= beta rho, and |phi(t)| <= e^{(Im t)^2/2}/sqrt(2 pi)).""" + br = beta * Iv.frac(rho) + return 1 + br * (br.sq() / 2).exp() / SQRT2PI + + +class Model: + """int_0^inf g(T) Phi(alpha - beta T) dT for g = sum const * poly(T) e^{-gamma T^2} [erf(kappa T)] [K_{al,g}(T)] (entire). On each + cell of width CELL the point Taylor coefficients of g at the centre are computed once and contracted with the cell moments, so + that an evaluation for given (alpha, beta) costs one product of coefficient vectors per cell. The remainder of the order-K + Taylor polynomial of the product is bounded by Cauchy's estimate at the centre, |c_n| <= M_g M_Phi / rho^n, summed geometrically; + beyond T_CUT the integrand is bounded by |g| (0 <= Phi <= 1, |erf| <= 1, |K_{al,g}| <= pi/(4 sqrt(al g))) and incomplete gamma + bounds. K_{al,g} is a certified cumulative layer (Math-#204 section 3): its value at a centre is the certified integral of + k = e^{-al w^2} E_g(w) over the full cells before plus the left half of the current cell, its higher Taylor coefficients are those + of k divided by n, and its disc bound is |K(s0)| + rho M_k.""" + + def __init__(self, terms, K=K_ORDER, cell=CELL, T=T_CUT, rho=RHO): + self.K, self.rho = K, rho + half = cell / 2 + ncell = int(T / cell) + if ncell * cell != T or Fr(half) >= Fr(rho): + raise ValueError("bad quadrature parameters") + mom = [Iv.frac(2 * half ** (n + 1) / (n + 1)) if n % 2 == 0 else Iv(ZERO) for n in range(K + 1)] + left = [Iv.frac((-1) ** n * half ** (n + 1) / (n + 1)) for n in range(K + 1)] # int_{-h}^0 x^n dx + rem_moment = Iv.frac(2 * half ** (K + 2) / (K + 2) / (1 - Fr(half) / Fr(rho))) + rho_pow = Iv.frac(Fr(rho) ** (K + 1)) + kinds = sorted({kk for *_, kk in terms if kk is not None}) + Kacc = {kk: Iv(ZERO) for kk in kinds} + self.Kvals = {kk: [] for kk in kinds} + self.cells, Rsum = [], Iv(ZERO) + for c in range(ncell): + s0 = (2 * c + 1) * half + Kco, Kbd = {}, {} + for al, g in kinds: + Ec = (PI / Iv.frac(g)).sqrt() / 2 + kc = [Ec * x for x in t_mul(gauss_coeffs(Iv.frac(al), s0, K), erf_coeffs(Iv.frac(g).sqrt(), s0, K))] + Mk = Ec * gauss_bound(Iv.frac(al), s0, rho) * erf_bound(Iv.frac(g).sqrt(), s0, rho) + Rk = (Mk / rho_pow * rem_moment / 2).hi # one-sided Cauchy remainder + if MUT == "remainder-dropped": + Rk = ZERO + lefti = sum((kc[n] * left[n] for n in range(K + 1)), Iv(ZERO)) + full = sum((kc[n] * mom[n] for n in range(0, K + 1, 2)), Iv(ZERO)) + Ks0 = Kacc[(al, g)] + (Iv(ZERO) if MUT == "cumulative-shift" else lefti) + Iv(Rk.copy_negate(), Rk) + Kacc[(al, g)] = Kacc[(al, g)] + full + Iv(CC.multiply(Decimal(-2), Rk), CC.multiply(Decimal(2), Rk)) + Kco[(al, g)] = [Ks0] + [kc[n] / (n + 1) for n in range(K)] + Kbd[(al, g)] = Iv(Ks0.mag()) + Iv.frac(rho) * Mk + self.Kvals[(al, g)].append((s0, Ks0)) + G, M = [Iv(ZERO)] * (K + 1), Iv(ZERO) + for const, poly, gam, kap, kk in terms: + pc = t_mul(poly_coeffs(poly, s0, K), gauss_coeffs(gam, s0, K)) + Mb = poly_bound(poly, s0, rho) * gauss_bound(gam, s0, rho) + if kap is not None: + pc = t_mul(pc, erf_coeffs(kap, s0, K)) + Mb = Mb * erf_bound(kap, s0, rho) + if kk is not None: + pc = t_mul(pc, Kco[kk]) + Mb = Mb * Kbd[kk] + G = [x + const * y for x, y in zip(G, pc)] + M = M + Iv(const.mag()) * Mb + Gm = [] + for j in range(K + 1): + acc = Iv(ZERO) + for i in range(j % 2, K + 1 - j, 2): + acc = acc + G[i] * mom[i + j] + Gm.append(acc) + self.cells.append((s0, Gm)) + Rsum = Rsum + M / rho_pow * rem_moment + self.Rsum = Rsum.hi + self.K_at_cut = Kacc + tail = Iv(ZERO) + for const, poly, gam, kap, kk in terms: # int_T^inf t^j e^{-g t^2} dt = (1/2) g^{-(j+1)/2} Gamma((j+1)/2, g T^2) + kfac = Iv(ONE) if kk is None else PI / (4 * (Iv.frac(kk[0]) * Iv.frac(kk[1])).sqrt()) + x = gam * Iv.frac(T * T) + for j, cj in enumerate(poly): + if cj == 0: + continue + a = Fr(j + 1, 2) + if a.denominator == 2: + xa1, gpow = ipow(x, int(a - Fr(1, 2))) * x.sqrt() / x, ipow(gam, int(a)) * gam.sqrt() + else: + xa1, gpow = ipow(x, int(a) - 1), ipow(gam, int(a)) + up = xa1 * (-x).exp() # Gamma(a, x) <= x^(a-1) e^-x (a <= 1), / (1 - (a-1)/x) (a > 1) + if a > 1: + up = up / (1 - Iv.frac(a - 1) / x) + tail = tail + Iv(const.mag()) * kfac * Iv.frac(abs(cj)) * up / gpow / 2 + self.tail = tail.hi + + def integral(self, alpha=None, beta=None): + tot = Iv(ZERO) + if alpha is None: + for _, Gm in self.cells: + tot = tot + Gm[0] + R = self.Rsum + else: + for s0, Gm in self.cells: + Wc = phi_aff_coeffs(alpha, beta, s0, self.K) + for j in range(self.K + 1): + tot = tot + Wc[j] * Gm[j] + R = (Iv(self.Rsum) * phi_aff_bound(beta, self.rho)).hi + if MUT == "remainder-dropped": + R = ZERO + R = CC.add(R, self.tail) + return tot + Iv(R.copy_negate(), R) + + +MODELS = {} + + +def model(m, gamma_out=None): + key = (m, gamma_out) + if key not in MODELS: + MODELS[key] = Model(g_terms(m, gamma_out)) + return MODELS[key] + + +Z_MEHTA = {2: 8 * SQRT2PI, 3: 48 * SQRT2 * PI, 4: 1536 * PI, 5: 46080 * PI * SQRTPI} # Mehta's integral (Math-#199 section 1; Z_4, Z_5: Math-#202, #204) + + +class Birth: + """I_B = kappa0/Z_m int_0^inf e^{-gamma_out T^2} F_m(T) [Phi(sl (b_+ - c T)) - Phi(sl (b_- - c T))] dT for the Gaussian pair + f ~ N(0, v), g = mu f + sigma Z (A = Q - g I): conditioning on f = b and integrating the GOE shift exactly (NOTE section 2). + Reference: v = 2/3, mu = 1, sigma = 0, i.e. kappa0 = sqrt(3/(m+3)), gamma_out = 3/(4m(m+3)), sl = sqrt((m+3)/2), c = 1/(m+3).""" + + def __init__(self, m, v=Fr(2, 3), mu=Fr(1), sig2=Fr(0)): + D = 1 + m * sig2 / 2 + A = 1 / (2 * v) + m * mu * mu / (4 * D) + Bc = mu / (2 * D) + g0 = Bc * Bc / (4 * A) + (Fr(0) if MUT == "tilt-gamma" else sig2 / (8 * D)) + self.m = m + self.gamma_out = Fr(1, 4 * m) - g0 + c = Bc / (2 * A) + if MUT == "trace-slope" and sig2 == 0: + c = Fr(1, m + 2) + self.sl = Iv.frac(2 * A).sqrt() + self.beta = self.sl * Iv.frac(c) + self.pref = 1 / Iv.frac(2 * v * A * D).sqrt() / Z_MEHTA[m] + self.M = model(m, None if self.gamma_out == F_pieces(m)[1] else self.gamma_out) + self.cache = {} + + def lam(self, b): + if isinstance(b, float): + return self.M.integral() if b > 0 else Iv(ZERO) + key = (b.lo, b.hi) + if key not in self.cache: + self.cache[key] = self.M.integral(self.sl * b, self.beta) + return self.cache[key] + + def I(self, B): + lo, hi = (x if isinstance(x, (float, Iv)) else Iv.frac(x) for x in B) + return self.pref * (self.lam(hi) - self.lam(lo)) + + +def gl_nodes(n): + xs, ws = [], [] + for i in range(1, n + 1): + x = math.cos(math.pi * (i - 0.25) / (n + 0.5)) + for _ in range(100): + p0, p1 = 1.0, x + for k in range(2, n + 1): + p0, p1 = p1, ((2 * k - 1) * x * p1 - (k - 1) * p0) / k + dp = n * (x * p1 - p0) / (x * x - 1) + dx = p1 / dp + x -= dx + if abs(dx) < 1e-15: + break + xs.append(x) + ws.append(2 / ((1 - x * x) * dp * dp)) + return xs, ws + + +def F5_direct_float(T, n=16): + """int_{T > w' > w > s > q > 0} P e^{-w'^2/80 - w^2/48 - s^2/24 - q^2/8} by nested Gauss-Legendre (float control of the exact + layers and of the cumulative K layer; P = prod mu_i^2 times the Vandermonde in the gaps).""" + xs, ws = gl_nodes(n) + tot = 0.0 + for x1, w1 in zip(xs, ws): + u = T * (x1 + 1) / 2 + a = (T - u) / 5 + for x2, w2 in zip(xs, ws): + w = u * (x2 + 1) / 2 + p1 = (u - w) / 4 + for x3, w3 in zip(xs, ws): + s_ = w * (x3 + 1) / 2 + p2 = (w - s_) / 3 + pre = w1 * w2 * w3 * (T / 2) * (u / 2) * (w / 2) * (s_ / 2) * math.exp(-u * u / 80 - w * w / 48 - s_ * s_ / 24) + for x4, w4 in zip(xs, ws): + q = s_ * (x4 + 1) / 2 + p3, p4 = (s_ - q) / 2, q + m1, m2, m3, m4, m5 = a, a + p1, a + p1 + p2, a + p1 + p2 + p3, a + p1 + p2 + p3 + p4 + gaps = (p1, p2, p3, p4) + vand = 1.0 + for i in range(4): + acc = 0.0 + for j in range(i, 4): + acc += gaps[j] + vand *= acc + tot += pre * w4 * (m1 * m2 * m3 * m4 * m5) ** 2 * vand * math.exp(-q * q / 8) + return tot + + +def I5_float(B, n=48, R=80.0): + """Float control of the outer T layer: Gauss-Legendre in T of e^{-3T^2/160} F_5(T) W_B(T) with F_5 from the float evaluation of + the exact pieces (math.erf and Simpson's K), independent of the certified Taylor/Cauchy quadrature and of the cumulative layer.""" + xs, ws = gl_nodes(n) + pieces = F_pieces(5)[0] + lo, hi = float(B[0]), float(B[1]) + + def Phi_f(x): + return 0.5 * math.erfc(-x / math.sqrt(2)) + tot = 0.0 + for k in range(4): # four panels of R/4 + a0 = k * R / 4 + for x, w in zip(xs, ws): + T = a0 + R / 8 * (x + 1) + W = Phi_f(2 * (hi - T / 8)) - Phi_f(2 * (lo - T / 8)) + tot += w * R / 8 * math.exp(-3 * T * T / 160) * F_float(pieces, T) * W + return math.sqrt(3 / 8) / (46080 * math.pi ** 1.5) * tot + + +# ---------------- the d = 6 jet of the reference kernel exp(-|z|^2/2), frame (u, w1, ..., w5) +def _e(*idx): + v = [0] * 6 + for i in idx: + v[i] += 1 + return tuple(v) + + +A_PAIRS = [(i, j) for i in range(1, 6) for j in range(i, 6)] # A11 A12 ... A55 (15 entries) +EVEN = [_e(), _e(0, 0)] + [_e(0, i) for i in range(1, 6)] + [_e(i, j) for i, j in A_PAIRS] # f, f_uu, f_uw1..5, A (22 entries) +ODD = [_e(0)] + [_e(i) for i in range(1, 6)] + [_e(0, 0, 0)] # f_u, f_w1..5, t = f_uuu +IDX_X, IDX_V = [0] + list(range(7, 22)), [1, 2, 3, 4, 5, 6] +DIAG_A = [1 + k for k, (i, j) in enumerate(A_PAIRS) if i == j] # positions of A11, ..., A55 in IDX_X +OFF_A = [1 + k for k, (i, j) in enumerate(A_PAIRS) if i != j] +N_VEC, DEGREE = 16, 10 # Lemma S: dimension of (f, A), degree of det(A)^2 + + +def dfact2(n): + r = 1 + while n > 1: + r *= n + n -= 2 + return r + + +def ref_cov(a, b): + """Cov(d^a f, d^b f) = (-1)^|a| d^(a+b) exp(-|z|^2/2) at 0 = (-1)^|a| prod_i (-1)^(al_i/2) (al_i - 1)!! (all al_i even).""" + al = [x + y for x, y in zip(a, b)] + if any(x % 2 for x in al): + return Fr(0) + v = Fr((-1) ** sum(a)) + for x in al: + v *= (-1) ** (x // 2) * dfact2(x - 1) + if MUT == "tau-cross" and sorted([a, b]) == sorted([_e(0, 0, 0), _e(0)]): + return Fr(0) + return v + + +def ref_matrix(idx): + return [[ref_cov(a, b) for b in idx] for a in idx] + + +def det_n(M): + """Determinant by Laplace expansion along the first row (n <= 6); Fractions or intervals.""" + n = len(M) + if n == 1: + return M[0][0] + if n == 2: + return M[0][0] * M[1][1] - M[0][1] * M[1][0] + tot = None + for j in range(n): + minor = [[M[i][k] for k in range(n) if k != j] for i in range(1, n)] + t = M[0][j] * det_n(minor) + tot = t if tot is None else (tot + t if j % 2 == 0 else tot - t) + return tot + + +def inv_n(M): + """Inverse by the adjugate (n <= 6); works for Fractions and for intervals.""" + n = len(M) + d = det_n(M) + cof = [[None] * n for _ in range(n)] + for i in range(n): + for j in range(n): + minor = [[M[r][c] for c in range(n) if c != j] for r in range(n) if r != i] + v = det_n(minor) + cof[i][j] = v if (i + j) % 2 == 0 else -v + return [[cof[j][i] / d for j in range(n)] for i in range(n)], d + + +def schur(C, X, V): + """Cov(X | V) = C_XX - C_XV C_VV^-1 C_VX for index lists into C (Fractions or intervals).""" + Vinv, detV = inv_n([[C[i][j] for j in V] for i in V]) + out = [] + for i in X: + row = [] + for j in X: + s = None + for p in range(len(V)): + for q in range(len(V)): + t = C[i][V[p]] * Vinv[p][q] * C[V[q]][j] + s = t if s is None else s + t + row.append(C[i][j] + s if MUT == "schur-sign" else C[i][j] - s) + out.append(row) + return out, detV + + +def sylvester_pd(M): + """Exact positive definiteness of a symmetric rational matrix (Gaussian elimination; all pivots positive).""" + n = len(M) + A = [list(map(Fr, r)) for r in M] + for k in range(n): + if A[k][k] <= 0: + return False + for i in range(k + 1, n): + f = A[i][k] / A[k][k] + for j in range(k, n): + A[i][j] -= f * A[k][j] + return True + + +C_EVEN_REF = ref_matrix(EVEN) +C_ODD_REF = ref_matrix(ODD) +CPRIME_REF, DETV_REF = schur(C_EVEN_REF, IDX_X, IDX_V) +LAMBDA0 = Fr(23, 100) + + +def shifted(M, lam): + return [[M[i][j] - (lam if i == j else 0) for j in range(len(M))] for i in range(len(M))] + + +# ---------------- image bound for the torus kernel in d = 6 +def image_bound6(L): + """E^(6)_L = 314496 (76 L^6 + 15) exp(-L^2/2): the shell |n|_inf = j carries (2j+1)^6 - (2j-1)^6 = 384 j^5 + 320 j^3 + 24 j <= 728 j^5 + points with |n|^2 <= 6 j^2, so sum_{n != 0} |n|^6 exp(-L^2 |n|^2/2) <= 157248 sum_j j^11 exp(-L^2 j^2/2) <= 314496 exp(-L^2/2) + (successive terms ratio at most 2^11 exp(-3L^2/2) < 1/2 for L >= 10); the contraction bound 76|x|^6 phi(x) (q <= 6, |x| >= 1), + |D^q phi(0)| <= 15 and the normalization S >= 1 are those of Math-#197 section 3 (dimension-free); the constant is Math-#204's E_6.""" + shells = 384 if MUT == "image-shells" else 728 + L = Fr(L) + return Iv.frac(2 * shells * 216 * (76 * L ** 6 + 15)) * (-Iv.frac(L * L / 2)).exp() + + +def half_power(x, twice): + r = ipow(x, twice // 2) + return r * x.sqrt() if twice % 2 else r + + +SANDWICH_VIOLATIONS = [] + + +def scaled(b, s): + if isinstance(b, float): + return b + b = b if isinstance(b, Iv) else Iv.frac(b) + return b * s if MUT == "window-scale" else b / s + + +def sandwich(birth, B, eps): + """Lemma S with n = 16, degree 10: (1-eps)^13/(1+eps)^8 I^ref_{B/sqrt(1-eps)} <= I'_B <= (1+eps)^13/(1-eps)^8 + I^ref_{B/sqrt(1+eps)}.""" + e = Iv(eps) if isinstance(eps, Decimal) else Iv.of(eps) + if e.hi >= 1: + raise ValueError("eps must be below 1") + up, dn = 1 + e, 1 - e + s_up, s_dn = up.sqrt(), dn.sqrt() + su, sd = (s_dn, s_up) if MUT == "sandwich-swap" else (s_up, s_dn) + lo, hi = B + upper = half_power(up, N_VEC + DEGREE) / half_power(dn, N_VEC) * birth.I((scaled(lo, su), scaled(hi, su))) + lower = half_power(dn, N_VEC + DEGREE) / half_power(up, N_VEC) * birth.I((scaled(lo, sd), scaled(hi, sd))) + if lower.lo > upper.hi: # impossible for a correct bracket; recorded and reported by rule SANDWICH_ORDERED + SANDWICH_VIOLATIONS.append(str(B)) + return Iv(upper.lo, lower.hi) + return Iv(lower.lo, upper.hi) + + +def scaled_exact(birth, B, eta): + """Exact I'_B for C' = (1 + eta) C'_ref: (1 + eta)^5 I^ref_{B/sqrt(1+eta)} (degree-10 homogeneity of det(A)^2).""" + s = (1 + Iv.frac(eta)).sqrt() + lo, hi = B + return ipow(1 + Iv.frac(eta), DEGREE // 2) * birth.I((scaled(lo, s), scaled(hi, s))) + + +def frob_eps(Cp, Cref): + """eps = ||C' - C'_ref||_F / lambda0, entrywise supremum of the deviation over each interval entry.""" + dev2 = Iv(ZERO) + n = len(Cref) + for i in range(n): + for j in range(n): + r = Iv.frac(Cref[i][j]) + c = Cp[i][j] if isinstance(Cp[i][j], Iv) else Iv.frac(Cp[i][j]) + d = max(abs(CC.subtract(c.hi, r.lo)), abs(CF.subtract(c.lo, r.hi))) + dev2 = dev2 + Iv(d).sq() + return (dev2.sqrt() / Iv.frac(LAMBDA0)).hi + + +class Space6: + """Every covariance entry of the 6-jet (other than Var f = 1, exact by the normalization) in [ref - E, ref + E]; the odd block + (G, t) is independent of the even block (f, V, A) for any even kernel ([LP] section 15).""" + + def __init__(self, E, birth): + E = Iv.of(E) + self.birth = birth + + def box(M): + n = len(M) + Bx = [[None] * n for _ in range(n)] + for i in range(n): + for j in range(i, n): + Bx[i][j] = Bx[j][i] = Iv.frac(M[i][j]) + Iv(E.hi.copy_negate(), E.hi) + return Bx + CE_ = box(C_EVEN_REF) + CE_[0][0] = Iv(ONE) + CO = box(C_ODD_REF) + G = [[CO[i][j] for j in range(6)] for i in range(6)] + Ginv, detG = inv_n(G) + c = [CO[6][j] for j in range(6)] + self.tau2 = CO[6][6] - sum((c[i] * Ginv[i][j] * c[j] for i in range(6) for j in range(6)), Iv(ZERO)) + self.pG = 1 / (8 * PI.sq() * PI * detG.sqrt()) # (2 pi)^(-3) det(Cov G)^(-1/2) + self.Cp, detV = schur(CE_, IDX_X, IDX_V) + self.pV = 1 / (8 * PI.sq() * PI * detV.sqrt()) + self.tau43 = root_iv(self.tau2.sq(), 3) + self.eps = frob_eps(self.Cp, CPRIME_REF) + + def k_integral(self, K): + lo, hi = K + + def s(k): + return k if isinstance(k, float) else 72 * Iv.frac(k * k) / self.tau2 + a = SHAPE() + return self.tau43 * C72 * (gamma_lower(s(hi), a) - gamma_lower(s(lo), a)) / (2 * SQRT2PI) + + def coefficient(self, B, K): + sphere = 8 * PI.sq() / 3 if MUT == "sphere" else PI.sq() * PI # |S^5| = pi^3 + return 144 * sphere * self.pG * self.pV * sandwich(self.birth, B, self.eps) * self.k_integral(K) + + +# ---------------- a non-scalar exact test of Lemma S: the tilted pair (f, g) +TILT = (Fr(17, 25), Fr(33, 34), Fr(2261, 86700)) # v = Var f, mu, sigma^2 with mu^2 v + sigma^2 = 2/3 + + +def tilted_cprime(v, mu, sig2): + """Covariance of (f, A11, A12, ..., A55) for A = Q - (mu f + sigma Z) I.""" + s = mu * mu * v + sig2 + C = [[Fr(0)] * 16 for _ in range(16)] + C[0][0] = v + for i in DIAG_A: + C[0][i] = C[i][0] = -mu * v + for j in DIAG_A: + C[i][j] = s + (2 if i == j else 0) + for i in OFF_A: + C[i][i] = Fr(1) + return C + + +# ---------------- Owen's T closed form of the Owen-type kinds, and the interval value of F_5 (rules OWEN_EXACT, INNER_LAYERS) +def K_owen(al, g, T): + """K_{al,g}(T) = int_0^T e^{-al w^2} E_g(w) dw = (pi/sqrt(al g)) [arctan(a)/(2 pi) - T_Owen(sqrt(2 al) T, a)], a = sqrt(g/al) > 1, + through the reflection T(h, a) = (Phi(h) + Phi(a h))/2 - Phi(h) Phi(a h) - T(a h, 1/a) (h >= 0) and arctan a = pi/2 - arctan(1/a); + Owen's series of the Math-#197 toolkit with N = 1400 terms (independent of the cumulative Taylor layer).""" + A = (Iv.frac(g) / Iv.frac(al)).sqrt() + h = (2 * Iv.frac(al)).sqrt() * Iv.frac(T) + ah = A * h + Th = (Phi(h) + Phi(ah)) / 2 - Phi(h) * Phi(ah) - owen_T(ah, 1 / A, N=1400) + return PI / (Iv.frac(al) * Iv.frac(g)).sqrt() * ((PI / 2 - atan_iv(1 / A, N=400)) / (2 * PI) - Th) + + +def F_iv(pieces, T): + """Interval value of F_m(T) from its closed-form pieces (the Owen-type kinds by K_owen); exact rational T.""" + tot = Iv(ZERO) + Ti = Iv.frac(T) + for key, poly in pieces.items(): + poly = trim(poly) + if not poly: + continue + pv = Iv.frac(sum(c * Fr(T) ** i for i, c in enumerate(poly))) + if key[0] == "1": + tot = tot + pv + elif key[0] == "e": + tot = tot + pv * (-(Iv.frac(key[1]) * Ti.sq())).exp() + elif key[0] == "E": + tot = tot + pv * (PI / Iv.frac(key[1])).sqrt() / 2 * erf_iv(Iv.frac(key[1]).sqrt() * Ti) + elif key[0] == "eE": + al, g = key[1], key[2] + tot = tot + pv * (-(Iv.frac(al) * Ti.sq())).exp() * (PI / Iv.frac(g)).sqrt() / 2 * erf_iv(Iv.frac(g).sqrt() * Ti) + else: + tot = tot + pv * K_owen(key[1], key[2], T) + return tot + + +# ---------------- main +S6, S21, S14 = Iv.frac(6).sqrt(), Iv.frac(21).sqrt(), Iv.frac(14).sqrt() +ATAN2 = PI / 2 - atan_iv(Iv.frac(Fr(1, 2))) +D3_CF = (50 * PI + 200 * ATAN2 - 228) / (9 * PI) # Math-#201 +D4_CF = (Iv.frac(Fr(6695, 54)) - Iv.frac(Fr(405, 32)) * S6 # Math-#202 + - (Iv.frac(Fr(1375, 72)) * S21 + Iv.frac(Fr(6695, 27)) * atan_iv(S21 / 7) + Iv.frac(Fr(405, 16)) * S6 * atan_iv(S14 / 14)) / PI) +R8 = [Fr(3, 1024) * x for x in (11243520, 0, -522240, 0, 28800, 0, -288, 0, 1)] +R6 = [Fr(-256, 729) * x for x in (405, 0, 405, 0, 45, 0, 1)] +STATED_F4 = {("e", Fr(1, 16)): [Fr(1, 31104) * x for x in (0, -991830528, 0, 18489408, 0, -346856, 0, 729)], # Math-#213, verbatim + ("E", Fr(1, 16)): R8, + ("eE", Fr(1, 48), Fr(1, 24)): R6, ("eE", Fr(1, 48), Fr(1, 6)): [2 * x for x in R6], + ("e", Fr(3, 16)): [Fr(512, 243) * x for x in (0, -297, 0, -42, 0, -1)]} +R10 = (7189875000, 0, 554203125, 0, -14962500, 0, 176750, 0, -700, 0, 1) +R9 = (0, -349281450000, 0, -34608500625, 0, 39341075, 0, -521235, 0, 729) +STATED_F5 = {("e", Fr(1, 5)): [Fr(256, 10546875) * x for x in (-131412150000, 0, -13096721250, 0, 7502375, 0, -174960, 0, 243)], + ("e", Fr(3, 40)): [Fr(-1, 337500000) * x for x in (-1076528332800000, 0, 241731011520000, 0, -113418584000, 0, + 4666529880, 0, 150240501)], + ("eE", Fr(1, 30), Fr(1, 24)): [Fr(256, 158203125) * x for x in R9], + ("eE", Fr(1, 30), Fr(1, 6)): [Fr(512, 158203125) * x for x in R9], + ("eE", Fr(1, 80), Fr(1, 16)): [Fr(-3, 500000000) * x for x in (0, -211803379200000, 0, 4362224640000, 0, -26148412800, 0, + 290979840, 0, 9372409)], + ("K", Fr(1, 80), Fr(1, 16)): [Fr(576, 9765625) * x for x in R10], + ("K", Fr(1, 30), Fr(1, 24)): [Fr(768, 9765625) * x for x in R10], + ("K", Fr(1, 30), Fr(1, 6)): [Fr(1536, 9765625) * x for x in R10]} +STAGES_F5 = [(Fr(1, 24), ["('1',)"]), + (Fr(1, 48), ["('1',)", "('E', Fraction(3, 16))", "('K', Fraction(1, 48), Fraction(1, 24))", + "('K', Fraction(1, 48), Fraction(1, 6))"]), + (Fr(1, 80), ["('1',)", "('E', Fraction(1, 5))", "('E', Fraction(3, 40))"])] # generated by the derivation and zero +OWEN_KINDS = ((Fr(1, 80), Fr(1, 16)), (Fr(1, 30), Fr(1, 24)), (Fr(1, 30), Fr(1, 6))) + + +def label(W): + def f(x): + return ("inf" if x > 0 else "-inf") if isinstance(x, float) else str(x) + return "[%s,%s]" % (f(W[0]), f(W[1])) + + +def within(iv, fv, tol): + f, t = Decimal(repr(fv)), Decimal(repr(tol)) + return CE.subtract(iv.lo, f) <= t and CE.subtract(f, iv.hi) <= t + + +def common_prefix(iv): + lo, hi = str(iv.lo), str(iv.hi) + i = 0 + while i < min(len(lo), len(hi)) and lo[i] == hi[i]: + i += 1 + return lo[:i] + + +def rel_width(iv): + return CC.divide(iv.width(), CC.multiply(Decimal(2), iv.lo)) + + +def main(): + global MUT, C_EVEN_REF, C_ODD_REF, CPRIME_REF, DETV_REF, LAMBDA0 + ap = argparse.ArgumentParser() + ap.add_argument("--mutant", choices=MUTANTS) + MUT = ap.parse_args().mutant + C_EVEN_REF, C_ODD_REF = ref_matrix(EVEN), ref_matrix(ODD) # rebuilt: they depend on the mutant + CPRIME_REF, DETV_REF = schur(C_EVEN_REF, IDX_X, IDX_V) + LAMBDA0 = Fr(1, 4) if MUT == "d5-floor" else Fr(23, 100) + t0 = time.time() + checks = {} + full, fullK = (-INF, INF), (Fr(0), INF) + band, bandK = (Fr(0), Fr(1)), (Fr(1, 2), Fr(2)) + + # the reference law, exactly + stated = [[Fr(0)] * 16 for _ in range(16)] + stated[0][0] = Fr(2, 3) + for i in DIAG_A: + stated[0][i] = stated[i][0] = Fr(-2, 3) + for j in DIAG_A: + stated[i][j] = Fr(8, 3) if i == j else Fr(2, 3) + for i in OFF_A: + stated[i][i] = Fr(1) + s2 = CPRIME_REF[0][0] + mean = [CPRIME_REF[i][0] / s2 for i in range(1, 16)] + condA = [[CPRIME_REF[i][j] - CPRIME_REF[i][0] * CPRIME_REF[0][j] / s2 for j in range(1, 16)] for i in range(1, 16)] + goe = [[(Fr(2) if i + 1 in DIAG_A else Fr(1)) if i == j else Fr(0) for j in range(15)] for i in range(15)] + Ginv, detG = inv_n([[C_ODD_REF[i][j] for j in range(6)] for i in range(6)]) + tau2_ref = C_ODD_REF[6][6] - sum(C_ODD_REF[6][i] * Ginv[i][j] * C_ODD_REF[6][j] for i in range(6) for j in range(6)) + checks["REFERENCE_EXACT"] = bool(CPRIME_REF == stated and mean == [-1 if i + 1 in DIAG_A else 0 for i in range(15)] + and condA == goe and DETV_REF == 3 and detG == 1 and tau2_ref == 6 + and tilted_cprime(Fr(2, 3), Fr(1), Fr(0)) == CPRIME_REF) + checks["LAMBDA_FLOOR"] = bool(LAMBDA0 == Fr(23, 100) and sylvester_pd(shifted(CPRIME_REF, LAMBDA0)) + and sylvester_pd(shifted(CPRIME_REF, Fr(231, 1000))) + and not sylvester_pd(shifted(CPRIME_REF, Fr(2312, 10000))) + and not sylvester_pd(shifted(CPRIME_REF, Fr(1, 4)))) + + # the exact inner layers: the generated cancellations, the stated closed form, m = 4 against Math-#213, and the interval value of + # F_5 (Owen's T for the three Owen-type kinds) against a direct float quadruple integral + pieces5, gamma5, stages5 = F_pieces(5) + nonzero = {k: trim(v) for k, v in pieces5.items() if trim(v)} + cancelled = sorted("%s" % (k,) for k, v in pieces5.items() if not trim(v)) + pieces4 = F_pieces(4)[0] + layer_ctrl = {} + for T in (3, 6, 10): + fi = F_iv(pieces5, Fr(T)) + layer_ctrl[str(T)] = {"F5_interval": fi.pair(), "F5_direct_float_n20": "%.15g" % F5_direct_float(float(T), 20)} + checks["INNER_LAYERS"] = bool(nonzero == STATED_F5 and stages5 == STAGES_F5 and gamma5 == Fr(3, 160) + and cancelled == STAGES_F5[-1][1] + and {k: trim(v) for k, v in pieces4.items() if trim(v)} == STATED_F4 + and all(within(Iv(Decimal(v["F5_interval"][0]), Decimal(v["F5_interval"][1])), + float(v["F5_direct_float_n20"]), 1e-10 * abs(float(v["F5_direct_float_n20"]))) + for v in layer_ctrl.values())) + + # the reference birth integral: the cumulative layer against Owen's T, D_5, coarse quadrature, m = 2, 3, 4, float control + ref5 = Birth(5) + owen = {} + ok = True + for kk in OWEN_KINDS: + rows = {} + for s0, v in [ref5.M.Kvals[kk][0], ref5.M.Kvals[kk][40]] + [(T_CUT, ref5.M.K_at_cut[kk])]: + w = K_owen(kk[0], kk[1], s0) + ok &= v.intersects(w) and v.width() < Decimal("1e-50") and w.width() < Decimal("1e-50") + rows[str(s0)] = {"cumulative_layer": v.pair(), "owen_T": w.pair()} + owen["K_%s_%s" % kk] = rows + checks["OWEN_EXACT"] = bool(ok) + I_full = ref5.I(full) + D5_IV = Iv(Decimal(D5_204[0]), Decimal(D5_204[1])) + checks["D5_CERTIFIED"] = bool(I_full.intersects(D5_IV) and I_full.width() < Decimal("1e-50")) + coarse = ref5.pref * Model(g_terms(5), K=8, cell=Fr(2), T=Fr(80), rho=Fr(6)).integral() + checks["QUAD_COARSE"] = bool(coarse.contains(D5_IV) and coarse.contains(I_full)) + ref2 = Birth(2) + cross = {} + ok = ref2.I(full).intersects(D2_EXACT) and ref2.I(full).width() < Decimal("1e-45") + for B in ((Fr(0), Fr(1)), (Fr(3), Fr(7, 2)), (-INF, Fr(0))): + x, y = ref2.I(B), I3_ref(*B) + ok &= x.intersects(y) and x.width() < Decimal("1e-45") + cross["m2 " + label(B)] = {"T_quadrature": x.pair(), "Math197_closed_form": y.pair()} + ref3 = Birth(3) + F4band = ref3.I(band) / D3_CF + ok &= ref3.I(full).intersects(D3_CF) and str(F4band.lo).startswith(F4B_209) and str(F4band.hi).startswith(F4B_209) + cross["m3 full"] = ref3.I(full).pair() + cross["m3 F^(4)_[0,1]"] = F4band.pair() + ref4 = Birth(4) + F5band = ref4.I(band) / D4_CF + ok &= ref4.I(full).intersects(D4_CF) and str(F5band.lo).startswith(F5B_213) and str(F5band.hi).startswith(F5B_213) + cross["m4 full"] = ref4.I(full).pair() + cross["m4 F^(5)_[0,1]"] = F5band.pair() + checks["CROSSCHECK_M2_M4"] = bool(ok) + fctrl = {} + ok = True + for B in ((Fr(0), Fr(1)), (Fr(1), Fr(2))): + fv, cv = I5_float(B), ref5.I(B) + ok &= within(cv, fv, 1e-8 * fv) + fctrl[label(B)] = {"float_eigenvalue_quadrature": "%.10g" % fv, "certified": cv.pair()} + checks["FLOAT_CONTROL"] = bool(ok) + + # image bound and epsilon + E10, E12, E24 = image_bound6(10), image_bound6(12), image_bound6(24) + checks["IMAGE_BOUND"] = bool(Decimal("4.60e-9") < E10.lo and E10.hi < Decimal("4.62e-9") + and Decimal("3.82e-110") < E24.lo and E24.hi < Decimal("3.83e-110") + and (Iv.frac(2048) * (-Iv.frac(150)).exp()).hi < HALF) + S0, S10, S12, S24 = Space6(0, ref5), Space6(E10, ref5), Space6(E12, ref5), Space6(E24, ref5) + checks["EPSILON"] = bool(S0.eps < Decimal("1e-55") and S10.eps < Decimal("1e-6") and S12.eps < Decimal("1e-15") + and S24.eps < Decimal("1e-55")) + + # Lemma S against exact scalar perturbations and an exact non-scalar (tilted) covariance + ok = True + for eta in (Fr(1, 50), Fr(-1, 50)): + for B in ((Fr(0), Fr(1)), (Fr(3), Fr(7, 2)), (Fr(-7, 2), Fr(-3))): + ok &= sandwich(ref5, B, Iv.frac(abs(eta))).contains(scaled_exact(ref5, B, eta)) + checks["SANDWICH_EXACT_SCALING"] = bool(ok) + tb = Birth(5, *TILT) + eps_t = frob_eps(tilted_cprime(*TILT), CPRIME_REF) + full_t = tb.I(full) + tilt = {"v_mu_sigma2": [str(x) for x in TILT], "eps": str(eps_t), "full_window": full_t.pair(), "windows": {}} + ok = full_t.intersects(I_full) and Decimal("0.108") < eps_t < Decimal("0.109") + for B in ((Fr(0), Fr(1)), (Fr(1), Fr(2))): + if not ok: # eps outside its certified range: no bracket is formed + break + ex, sw, rf = tb.I(B), sandwich(ref5, B, eps_t), ref5.I(B) + ok &= sw.contains(ex) and not ex.intersects(rf) + tilt["windows"][label(B)] = {"tilted_exact": ex.pair(), "sandwich": sw.pair(), "reference": rf.pair()} + checks["TILTED_EXACT"] = bool(ok) + + # the coefficient: reference, every L >= 10, L = 24 + C6_REF = GAMMA76 * CBRT32 * I_full / (64 * SQRT3 * PI.sq() * PI * SQRTPI) # c_{6,ref} = Gamma(7/6) (3/2)^(1/3) D_5/(64 sqrt3 pi^(7/2)) + b_windows = [(Fr(0), Fr(1)), (Fr(-1), Fr(1)), (Fr(-2), Fr(2)), (Fr(0), INF), (-INF, Fr(0)), full] + table = {} + for B in b_windows: + F6 = ref5.I(B) / I_full + for K in (bandK, fullK): + GK = G_ref(K) + table["B=%s K=%s" % (label(B), label(K))] = { + "F6_B": F6.pair(), "G_K": GK.pair(), "c6_ref": (C6_REF * F6 * GK).pair(), + "c6_torus_L_ge_10": S10.coefficient(B, K).pair(), "c6_torus_L_24": S24.coefficient(B, K).pair()} + F6_band = ref5.I(band) / I_full + c_band_ref = C6_REF * F6_band * G_ref(bandK) + c_band_0, c_band_10, c_band_24 = S0.coefficient(band, bandK), S10.coefficient(band, bandK), S24.coefficient(band, bandK) + c_full_0, c_full_10, c_full_24 = S0.coefficient(full, fullK), S10.coefficient(full, fullK), S24.coefficient(full, fullK) + checks["FACTORIZATION"] = bool(c_band_0.intersects(c_band_ref) and c_full_0.intersects(C6_REF)) + checks["C6_CONSISTENT"] = bool(str(C6_REF.lo).startswith(C6_REF_204) and str(C6_REF.hi).startswith(C6_REF_204) + and c_full_24.intersects(C6_REF) and c_full_10.contains(C6_REF)) + checks["NESTING"] = bool(c_band_10.contains(c_band_ref) and c_band_10.intersects(c_band_24) + and c_band_10.width() > c_band_24.width() * Decimal(10) ** 30 + and all(Iv(Decimal(v["c6_torus_L_ge_10"][0]), Decimal(v["c6_torus_L_ge_10"][1])).contains( + Iv(Decimal(v["c6_ref"][0]), Decimal(v["c6_ref"][1]))) for v in table.values())) + rel10 = rel_width(c_band_10) + worst = max(rel_width(Iv(Decimal(v["c6_torus_L_ge_10"][0]), Decimal(v["c6_torus_L_ge_10"][1]))) for v in table.values()) + checks["WIDTHS"] = bool(worst < Decimal("2e-5") and c_band_24.width() < c_band_24.lo * Decimal("1e-40")) + checks["SANDWICH_ORDERED"] = bool(not SANDWICH_VIOLATIONS) + checks["LIBRARY_EXACT"] = bool(library_exact()) + pinned = {"c_6_ref": C6_REF, "D_5": I_full, "F6|[0,1]": F6_band, "c6|[0,1]x[1/2,2] ref": c_band_ref, "pi": PI} + checks["PINNED"] = bool(all(str(pinned[k].lo).startswith(p) and str(pinned[k].hi).startswith(p) for k, p in PINNED.items())) + if len(checks) != 19: + raise ValueError("rule count changed") + passed = all(checks.values()) + + out = { + "object": "CL-C8-WINDOW-COEFFICIENT-D6-20261001-v1", + "scientific_effect": "NONE", + "certified": True, + "mutant": MUT, + "passed": passed, + "checks": checks, + "c8_band": {"B": "[0,1]", "K": "[1/2,2]", + "c6_reference": c_band_ref.pair(), "c6_reference_digits": common_prefix(c_band_ref), + "F6_B": F6_band.pair(), "c6_torus_every_L_ge_10_every_frame": c_band_10.pair(), + "c6_torus_L_24": c_band_24.pair(), "relative_half_width_L_ge_10": str(rel10)}, + "full_window": {"c_6_ref": C6_REF.pair(), "c6_torus_every_L_ge_10": c_full_10.pair(), "c6_torus_L_24": c_full_24.pair(), + "Math204_c_6_ref_digits": C6_REF_204}, + "windows": table, + "worst_relative_half_width_L_ge_10": str(worst), + "transfer": {"E6_10": E10.pair(), "E6_12": E12.pair(), "E6_24": E24.pair(), + "eps_L_ge_10": str(S10.eps), "eps_L_ge_12": str(S12.eps), "eps_L_24": str(S24.eps), "eps_E0": str(S0.eps), + "lambda_floor": str(LAMBDA0), "lambda_min_exact": "(6 - sqrt(92/3))/2", + "tau2_L_ge_10": S10.tau2.pair(), "pG_pV_L_ge_10": (S10.pG * S10.pV).pair()}, + "reference": {"C_prime": [[str(x) for x in r] for r in CPRIME_REF], "A_given_f_mean": [str(x) for x in mean], + "tau2": str(tau2_ref), "det_Cov_V": str(DETV_REF), + "D_5_quadrature": I_full.pair(), "D_5_Math204": list(D5_204), "coarse_quadrature": coarse.pair(), + "F_5_pieces": {"%s" % (k,): [str(x) for x in v] for k, v in sorted(nonzero.items(), key=lambda kv: str(kv[0]))}, + "layer_cancellations": [[str(c), z] for c, z in stages5], + "quadrature_remainder_sum": str(ref5.M.Rsum), "quadrature_tail": str(ref5.M.tail)}, + "owen_type_layer": owen, + "layers_interval_vs_float": layer_ctrl, + "crosschecks": cross, + "float_control": fctrl, + "tilted": tilt, + "constants": {"pi": PI.pair(), "Gamma(7/6)": GAMMA76.pair(), "(3/2)^(1/3)": CBRT32.pair()}, + "exact_relations": { + "window coefficient": "c^(6)_{B,K} = 144 int_{S^5} p_G(0) p_V(0) I_B(u) J_K(u) dsigma(u), I_B = E[1{f in B} det(A)^2 1{A<0} | V=0] ([LP] (11.3), section 15)", + "reference": "(f, A) | V=0: f ~ N(0, 2/3), A = Q - f I, Q the 5x5 GOE; I^ref_B = (sqrt(3/8)/(46080 pi^(3/2))) int_0^inf e^{-3T^2/160} F_5(T) [Phi(2(b_+ - T/8)) - Phi(2(b_- - T/8))] dT", + "Owen-type kinds": "K_{al,g}(T) = int_0^T e^{-al w^2} E_g(w) dw = (pi/sqrt(al g)) [arctan(a)/(2 pi) - T_Owen(sqrt(2 al) T, a)], a = sqrt(g/al)", + "sandwich": "(1-eps) C'_ref <= C' <= (1+eps) C'_ref => (1-eps)^13/(1+eps)^8 I^ref_{B/sqrt(1-eps)} <= I'_B <= (1+eps)^13/(1-eps)^8 I^ref_{B/sqrt(1+eps)}", + "epsilon": "eps = ||C' - C'_ref||_F / lambda0, lambda0 = 23/100 <= lambda_min(C'_ref) = (6 - sqrt(92/3))/2", + "image bound": "E^(6)_L = 314496 (76 L^6 + 15) exp(-L^2/2)", + "c_6_ref": "Gamma(7/6) (3/2)^(1/3) D_5 / (64 sqrt3 pi^(7/2)) (Math-#204)", + }, + } + sys.stdout.write(json.dumps(out, indent=1) + "\n") + sys.stderr.write("seconds %.1f\n" % (time.time() - t0)) + sys.exit(0 if passed else 1) + + +if __name__ == "__main__": + main()