tail.py
Checks: READING: The carry bound rule.
Cited on: READING.
To run it: save tail.py, with edge.py, every.py, heat.py, majorant.py, real_support.py, remainder.py, uniform.py beside it; needs mpmath, numpy, scipy.
"""
tail.py -- the real support's margin at every P >= P* = 399 2^25 10^20 =
1.3388e30, by one closed form in P.
QUESTION. box.py reads the margin at d = P + 1 = 0 mod 4 on boxes of P
from 399 to P*, where its range stops. Is there a closed form R(P) >= P
R_maj(P) at every P >= P*, with budget(P) - F(c(P), P) - X - R(P)/P > 0
at every such P? Notation as box.py's head (D1..D10), majorant.py (W1,
W2), uniform.py (X1, X2).
THE ANSWER. At every P >= P*,
P R_maj(P) <= C_K + L/P^2 + 8e-40 T(P) + W(P) + P past_T,
C_K = 7809.1, L = 1.30e8 (C6), each piece over P at most its value at
P* (THE MARGIN), so the margin is least at P*, 0.2497859, and positive
at every P >= P*. Conditions: W1 and W2 holding at P = P0 (majorant.py's
hand derivations), the pieces K1..K7 by hand as written below, and X <=
2.06e-4 past 1e6 by hand (variation.py). C_K and L are sums in mpmath at
30 digits, not interval arithmetic.
THE SCALE. With F <= 1e-5 (K1), X <= 2.06e-4 and budget > 1/4, the
margin is positive while R_maj(P) < 0.2497, i.e. P R_maj < 3.3e29 at P*.
A piece of P R_maj growing no faster than polylog(P) needs only a crude
constant: C_K <= 1e20 costs the margin 7.5e-11. Only pieces linear in P
(or faster) are owed sharply, and one is: P past_T. Every other
P-dependence is a floor made small by K(P), or the tau range [1/c(P),
T(P)] widening at both ends, paid by integrals that converge over (0,
inf).
THE PIECES (P-dependence named and signed)
(K1) c(P) := ceil(16 ln(1600 P)) + 1 (1229 at P*). By hand: g >=
(c/16)(1 - eps), g' >= (c/(4P))(1 - eps), eps < c/P, and 1/(1 -
e^-x) <= 1 + 1/x, so F <= 4 e^-g (1 + 1/g') <= (1/400)(1/P +
4/c)(1 + delta) <= 8.2e-6 at P*, falling in P (F_box reads
7.62e-6 at P*, 5.79e-6 at 1e40, 2.35e-6 at 1e100). DERIVED.
(K2) T(P) = 1 + ln(P/399)/(2 pi^2): past_T(P, T(P)) -> L_T = 8 399^2
e^(-4 pi^2)/(4 pi^2) = 2.30896e-13 (the k = 1 root: e^(-4 pi^2 T)
= e^(-4 pi^2)(399/P)^2 against 1/(1 - e^lr) ~ P^3/(4 pi^2)). Read:
uniform.past_T/L_T = 1.368 at 399, 1.0142 at 1e4, 1 + 2.2e-8 at
1e10, 1 +- 1.3e-13 at 1.3388e26, 1.3388e30, 1e40. In R_maj it is
past_T itself: 2.31e-13, P-free. By hand: d = P + 1, phi_k = 2 pi
k/d, x_k = ln Q_k = ln(1 + min(8 sin^2(phi_k/2), 2)), beta = T
P^3/(2d), v_k = x_k/(2d); the Q-bound summed over m >= T P^3 is
(8/P) sum_(k=1..d/2) e^(-beta x_k)/(1 - e^(-v_k)) (uniform.past_T
reads it, larger past k = 2000 by lq_beyond).
k = 1, z = pi/d: sin z >= z(1 - z^2/6) and ln(1 + y) >= y - y^2/2
give x_1 >= (8 pi^2/d^2)(1 - eta), eta = 13 pi^2/(3 d^2), so beta
x_1 >= 4 pi^2 T (P/d)^3 (1 - eta) >= 4 pi^2 T(1 - 3/d - eta), and
e^(-4 pi^2 T) = e^(-4 pi^2)(399/P)^2; v_1 in [(4 pi^2/d^3)(1 - eta),
4 pi^2/d^3] and 1 - e^-v >= v(1 - v/2). So the k = 1 term is at
most L_T (d/P)^3 e^(4 pi^2 T(3/d + eta))/((1 - eta)(1 - 2 pi^2/d^3)).
k >= 2, phi_k <= pi/3: sin t >= 3t/pi on [0, pi/6] and ln(1 + y) >=
y ln3/2 on [0, 2] give x_k >= (9 ln3/pi^2) phi_k^2 >= 4 pi^2 k^2/d^2,
so beta x_k >= a k^2, a = 2 pi^2 T (P/d)^3, and 1/(1 - e^-v) <= 1 +
1/v <= 1 + d^3/(2 pi^2 k^2); with k^2 >= 4 + 4(k - 2) they add at
most (8/P) e^(-4a)(1 + d^3/(8 pi^2))/(1 - e^(-4a)), over L_T eps_2
= ((d/P)^3/2 + 4 pi^2/P^3) e^(-4 pi^2 T(1 - 6/d))/(1 - e^(-4a)).
phi_k > pi/3: x_k >= ln 3, at most d/2 terms, eps_3 = (d/2)(1 +
2d/ln3)(4 pi^2/P^3) e^(4 pi^2 T - T P^3 ln3/(2d)).
uniform.past_T's own terms past k = 2000 (R_maj's): x read as lb =
min(1.196 phi_2001^2, 0.846), at most d/2 of them, eps_4 = (d/2)(1 +
2d/lb)(4 pi^2/P^3) e^(4 pi^2 T - beta lb), beta lb >= 4.2e7 T past P
= 1e4 (0.423 T P^3/d while 0.846 binds, 9.45e7 T (P/d)^3 past d =
1.5e4).
So uniform.past_T(P, T(P)) <= L_T(1 + eps_T(P)), 1 + eps_T = (k = 1
factor) + eps_2 + eps_3 + eps_4, each falling in P past 399 (T/d, T
eta and e^(-4 pi^2 T) fall, beta lb rises on both branches of its
min): eps_T(P*) about 3.7e-28.
In R_maj: past_T <= 2.30897e-13 (1 + eps_T(P*)) at every P >= P*.
(K3) The floors, every P-dependent term of B: the walk P^3 e^-g/(1 -
e^(-P g')), the N sum P(...)e^(-K^2/(4 tau)), phi's modes past K
P sum 4 pi k e^(-4 pi^2 k^2 tau), the root tail 4P^2(sum_(k>K)
e^(-8.72 tau k^2) + (P/4 + 1)e^(-0.546 tau P^2)). Images below
1/4 with K_img(P) = ceil(sqrt(4 ln P + 100)) + 2, each floor
there <= 1e-40/P; modes from 1/64 (K6 reads images to 1/16, modes
from 1/64), each floor <= 1e-40/P below 1/4 and, above it, falling
as K6 (iv)'s lambda K_lo >= 12.5 e^(-4 pi^2 tau) does (lambda K_lo
at T(P*) is 8.0e-72 against 1e-40/P = 7.5e-71, below it past tau
4.15, so a flat floor would not hold o = 0 alone): each floor <= 1e-40
e^(-4 pi^2 tau)/P^2.
K_mode(P) = ceil(sqrt(7.34(3 ln P + 95))) (by hand, worst at tau =
1/64), 48 at P*. Counted at EVERY class, the floors add at most
8e-40 T(P) to P R_maj. Bad needs P|K_x| <= b + fl <= 2 max(b, fl),
and P|K_x| <= 2 fl holds at o = 0 alone (K6 (iv)), where P|K_x| <=
2b holds too, so [Bad](b + fl) <= b[P|K_x| <= 2b] + fl.
By hand, tau < 1/4, K = K_img(P), J = KP, n0 = 2 tau P^2 < P^2/2:
g = J^2/(2(n0 + BJ)) >= K^2/(1 + 2BK/P) >= K^2 - 2.01 K^3/P, and K^2
>= 4 ln P + 100, so P^3 e^-g <= e^(-100 + 2.01 K^3/P)/P; g' >= J/(2(n0
+ BJ)), so P g' >= K/(4 tau + 2BK/P) >= 0.99 K: the walk <= 3.8e-44/P.
The N sum, P(4 pi tau)^(-1/2) sum_(N >= K)(N + 1) tau^-1 e^(-N^2/(4
tau)), rises in tau while K^2/(4 tau) > 3/2: at tau = 1/4 it is at
most 2.26 P (K + 1) e^(-K^2)(1 + e^(-2K)) <= 2.3 (K + 1) e^-100/P^3.
tau >= 1/64, K = K_mode(P), 8.72 K^2/64 >= 1.00007(3 ln P + 95):
the root tail's first term 4P^2 e^(-8.72 tau (K + 1)^2)/(1 -
e^(-17.4 tau(K + 1))) falls in tau, at 1/64 at most 4.01 e^-95/P =
2.2e-41/P; over e^(-4 pi^2 tau)/P^2 it is at most 4.1 P^4 e^(-(8.72(K
+ 1)^2 - 4 pi^2)tau), largest at tau = 1/4, there below e^(-44 ln P
- 1400); its second, 4P^2(P/4 + 1) e^(-0.546 tau P^2), and phi's modes
past K, at most P 4 pi (K + 1) e^(-4 pi^2 tau (K + 1)^2)(1.01) with
4 pi^2 (K + 1)^2/64 >= 13.5 ln P + 430, lie far below. So four
floors, fl <= 4e-40/P (tau < 1/4), 4e-40 e^(-4 pi^2 tau)/P^2 (tau >=
1/4), and over the T P^3 values of m at weight P/P^3 each, 2 fl adds
at most 8e-40 T(P) to P R_maj.
(K4) The images' arm, P-free: b_img(x, tau) = sum_(|k| <= K_img)
G(x + k, tau) Abar(|x + k|/(2 tau), tau) (K6 sums every k, above
it), A read at P0 = P (W1 holds at P = P0), Abar = 1.75(1.5 a + a^3
tau) + 1e-20. By hand: core images have a <= (K + 1)c/2, a/P <=
1.1e-26 at P*, so dbar, ybar <= 1.01 a/P; |l| <= 1.5a + a^3 tau;
dropping Hbar's subtracted terms, Hbar <= 1.25 a^2 (1 + e), and
|Lbar| <= (7/6) a^4 tau (1 + e) + (1 + e)/(6 tau), e <= 1e-25; and
1.25 a^2 <= (a/P)(1.5a) P, (7/6) a^4 tau = (7/6)(a/P) a^3 tau P, so
U, |V| <= (1.5a + a^3 tau)(1 + 2/P + 1.17 a/P) + 1.01/(6 tau P).
U <= (K + 1)^3 c^2/4 <= P (4.6e9 at P*), so f(U) <= (e - 1)U, and
A <= max(f(U), |V|) <= 1.7183(1 + 1e-25)(1.5a + a^3 tau) +
1.75 c(P)/(6P): below Abar, c/(6P) 1.5e-28 at P*. Each P-term
(a/P, c/P, (K + 1)^3 c^2/P) falls in P. Core: K_img + 1 <=
P/(2c(P)). DERIVED.
(K5) The modes' arm, P-free: b_mode = sum_(k >= 1) env_k Abar_k(tau)
bounds sum_(k <= K_mode(P)) at every P >= P*. W2 read at P0 = P
(it holds at P = P0), theta = 2 pi k/P <= 1e-6 for every k <=
K_mode(P) at every P >= P* (W2 asks theta0 <=
0.35, and variation.py's Q6 the disc's root). By hand at theta <=
1e-6: Ghat(r) <= 1.0834 r^4 (W2's exact |g_4| = 13/12, Cauchy past
it, 30.5 r^5); sigma = 1.01 is a fixed point's radius (left side
<= 1 + 1.0001e-6); mubar = 1.01 theta^2/P; rhobar <= 2.03
theta^3/P + 1.0835 theta^4, zeta0 = rhobar tau P^2 <= 1e-6; E_v <=
theta^3/6 + 1.02 theta^2/P; |v|, q <= (1 + 2e-6)theta; |E_D'| <=
1.50001 theta^2 + 2.02 theta/P; |e| <= 2.1667 theta^2 + 3.05
theta/P. So E_c = P|e| is O(k^2)/P, P(e^zeta0 - 1) O(k^4 tau)/P, and
A_k <= Abar_k := |beta1 - beta2| + C_k(tau)/P*,
C_k = sqrt(beta1 beta2^3/3) + 2.167(2 pi k)^2 + 3.06(2 pi k)
+ (2.04(2 pi k)^3 + 1.084(2 pi k)^4) tau,
each 1/P term falling in P; at beta1 = beta2 (tau = 3/(8 pi^2 k^2))
C_k is the whole of it, about 201.03 k^2 + 38.5 k against W2's own
201 k^2 + 37.7 k. DERIVED.
(K6) The window, one inequality on the whole period, no zones. K(x,
tau) = theta_3(pi x, q), q = e^(-4 pi^2 tau), and Jacobi's triple
product (known) gives -K_x = 4 pi K sin(2 pi x) sum_(m >= 1)
q^(2m-1)/|1 + q^(2m-1) e^(2 pi i x)|^2, each denominator at most (1
+ q)^2, so at every x
|K_x(x, tau)| >= lambda(tau) |sin 2 pi x| K(x, tau),
lambda = 4 pi q/((1 + q)^3 (1 - q)) (about 1/(8 pi tau) at 0).
A class x = x0 + u, u = o/P, |u| <= 1/4, has |sin 2 pi x| >= 4|u|,
so P|K_x| >= 4|o| lambda K(x). With K3's split, Bad needs P|K_x| <=
2b(x), so |o| <= N := r/(2 lambda), r = sup_x b/K, P-free.
(i) Images, tau <= 1/16, b = sum_k G0(x + k) w(a_k), G0 the heat
Gaussian, a_k = |x + k|/(2 tau), w(a) = a Abar(a, tau), K4's Abar
needing only: w rising, w(ca) <= c^4 w(a) at c >= 1. With y the
nearest image (|y| <= 1/2), K >= G0(y); image y + j, j >= 1 (j <=
-1 alike), z = y + j in [j - 1/2, j + 1/2], weighs w(z/(2 tau))
e^(-j(2z - j)/(4 tau)) <= (2z)^4 w(1/(4 tau)) e^(-j(2z - j)/(4
tau)), falling in z past 8 tau/j <= 1/2: so r <= rho w(1/(4 tau)),
rho = 1 + 2 sum_(j>=1) (2j - 1)^4 e^(-j(j - 1)/(4 tau)) (3 + 2e-12
at 1/64). The truth at x = 1/2 is about w(1/(4 tau)): r is 3 high,
4 lambda about pi/(8 tau) low against K_xx(1/2)/K(1/2) (480 at
1/c(P*)), so N is about 0.26/tau^2 against the pair's 8 tau Abar:
crude, and C_K affords it.
(ii) Modes, tau >= 1/64, b_mode x-free (K5): N = b_mode/(2 lambda
K_lo), K(x) >= K(1/2) >= K_lo = max(2 G0(1/2), 1 - 2q).
(iii) Bad classes lie in |u| <= u_w := (N + 1/2)/P_eff, P_eff =
max(P*, P_min(tau)), P_min(tau) = e^(1/(16 tau) - 1/8)/1600 the
least P whose [1/c(P), T(P)] can hold tau (c <= 16 ln(1600 P) + 2);
u_w/tau <= 4e-22 at every tau (3.538e-22 at most, at 1/1229.95,
where P_min = P*). There, with z
over |x0 + Z| less 0, two of each,
bloc_x0 := sup b <= (1 + 1e-20) (4 pi tau)^(-1/2) ([x0 = 0]
w(u_w/(2 tau))
+ 2 w(1/(4 tau)) sum_z (2z)^4 e^(-z^2/(4 tau))),
at x0 = 1/2 the pair's b(1/2) = 2 G(1/2) Abar(1/(4 tau)) itself.
(iv) The o = 0 clause (K3): |o| >= 1/2 gives P|K_x| >= 2 lambda
K_lo, so o = 0 alone has P|K_x| <= 2 fl once lambda K_lo > fl: 4e-7/P
against 4e-40/P below 1/4 (e^(-1/(16 tau)) >= e^(-1/8)/(1600 P) at
tau >= 1/c(P); lambda falls in tau, d ln/dq = (1 - 2q + 3q^2)/(q(1
- q^2)) > 0, 6.4e-4 at 1/4), 12.5 e^(-4 pi^2 tau) against 4e-40
e^(-4 pi^2 tau)/P^2 above. So no mv, m_far, L1 or far-zone cut is
owed. Counted with the wrong-signed class (its 2 UB; its 2P|K_x|
by K7's W(P)), at most 2N + 2 classes a
side: S_img = 2(2N + 2)(bloc_0 + bloc_1/2), S_mode = 4(2N + 2)
b_mode, the smaller on [1/64, 1/16]. Monotone in P: only u_w
moves, falling; N, lambda, K_lo are P-free.
(K7) The q-lattice: P R_maj = sum_q P^-2 (the class sum at tau_q, P^3
units), tau_q of step P^-2 in [1/c(P), T(P)), and S read at
P_eff(tau) bounds the class sum at every P whose range holds tau
(only u_w moves, falling in P). Cells [t1, t2) gapless on [tau0,
tau1] = [1/2000, 40], cut at 1/64 and 1/16, each holding at most
(t2 - t1)P^2 + 1 points: sum <= (t2 - t1 + P^-2) S_sup, S_sup each
monotone factor at its worse end: w(1/(4 tau)), (4 pi tau)^(-1/2)
and q^(k^2) falling, rho, 1/lambda, e^(-z^2/(4 tau)) and C_k
rising, w(a, tau) rising in both, a_w = (N + 1/2)/(2 tau P_eff)
at (N_sup, t1, P_eff(t2)), |beta1 - beta2| convex (larger end),
G0(1/2) unimodal (max at 1/8: lesser end), 1 - 2q rising; sums cut
at z, j <= 20 and k <= 60 (each remainder under e^-1665 its
first term; S carries them at most squared, so C_K and L rise by
under a factor 1 + e^-1600). Tails,
by hand: tau <= tau0, w(1/(4 tau)) <= 0.0072 tau^-3, 1/lambda <=
25.65 tau (q >= 0.98, 1 - q <= 4 pi^2 tau), rho <= 3.0001, so N <=
0.278 tau^-2, 2N + 2 <= 0.557 tau^-2; bloc_1/2 <= 0.0040627
tau^-3.5 e^(-1/(16 tau)), bloc_0 <= 1e-18 tau^-3.5 e^(-1/(16
tau)) (u_w <= 505 tau^-2
e^(-1/(16 tau)) by P_min, w(a) <= 1.0001e-20 a there): S <= f0 =
0.0046 tau^-5.5 e^(-1/(16 tau)), rising below 1/88, int_0^tau0 f0
= 0.0046 16^4.5 Gamma(4.5, 125). tau >= tau1: Abar_1 <= 40
tau^1.5 (Abar_1 below it from 38.4), the modes past 1 under
q^3, lambda >= 4 pi q(1 - 3q), K_lo >= 1 - 2q, so S <= 16 pi q
(Abar_1^2 + 2 Abar_1)(1 + e^-1000), 16 pi 1600 = 80424.8, the 2
Abar_1 inside the slack 75.2 tau^3 from tau 14.2: S <= f_inf = 80500
tau^3 e^(-4 pi^2 tau), falling.
Each tail's lattice sum is its integral plus P^-2 f at its end.
The wrong-signed class (odd P): 2P|K_x(1/2 - 1/(2P))| <= sup_x
|K_xx| <= F(tau) = sum_k 8 pi^2 k^2 q^(k^2), falling, so its sum
is at most int F from 1/c - P^-2 = 2 sum_k q^(k^2) there <= W(P) =
1/(2 sqrt(pi (1/c(P) - P^-2))), 9.9 at P*: polylog, in place of
the limit's int Gamma <= 1. So
P R_maj(P) <= C_K + P^-2 L + 8e-40 T(P) + W(P) + P past_T,
C_K = sum_cells (t2 - t1) S_sup + both tails' integrals, L =
sum_cells S_sup + f0(tau0) + f_inf(tau1), P-free.
The limit's first orders put int S near 18 (uniform.py C_R); S
itself peaks near 1/88 (f0's shape), so C_K sits decades above.
THE MARGIN. Each piece of the closed form over P is at most its value
at P* (T grows like ln P, W like sqrt(ln P); c(P) steps by 1, first at
1.0029 P*, raising W by 1.0004, its steps e^(1/16) apart in P after; at
each step F_K1's 4/(c(1 - c/P)) falls by about c/(c + 1), outweighing
its e^(c^2/(16P))'s rise, under 1e-27), so the margin is least at P*:
budget - 8.2e-6 - 2.06e-4 - C_K/P* - 2.31e-13 > 0.2497 for any C_K <
1e20.
THE CONTROLS, each frozen before its code. The closed form must
dominate the TRUE class sum, not box.py's U/Pa, which carries its box's
ratio (Pb past_T/Pa = 1e4 L_T on the last box).
(C1) S(tau) against the true class sum at P = P* (uniform.py's
summand rewritten per (x0, o): u = o/P carried apart from x0 and
K_x(x0 + u) evaluated from u (the Fourier form's sin(2 pi k(1/2 +
u)) = (-1)^k sin(2 pi k u); below tau = 0.05 the cancelling image
pair combined analytically, or extended precision), and the sign s
and img_B's x and j formed from (x0, o) too: in doubles an absolute
residue near 6.7e29 has spacing 1.4e14, and even 0.5 + 400/P* reads
0.5, so K_x would read 0 and every class Bad; the classes within 400
of each sign change, majorants at P0 = P*) at the centres of every
third of 300 geometric tau cells on [1/c(P*), T(P*)]: largest
ratio truth/S at most 1, and no Bad class at 300 <= |o| <= 400
(the scan wide enough). K6 has no far zone, so C1 sees it whole.
FROZEN BEFORE CODE: the truth is the closed form's own B, images
|k| <= K_img(P*) = 22 at W1's A (P0 = P*, the k = 0 image's xi
from u, the rest from x0: B moves by under 1e-24 relative), modes
k <= K_mode(P*) = 48 at W2's A_k (P0 = P*), the smaller where both
run; Bad iff s phi <= B + fl (K3's fl), each Bad class 2(b + fl);
the wrong-signed class's 2P|K_x| (odd P alone) printed apart over
Gamma (K7's W(P), not S). P = P* (even, every class right-
signed or s phi = 0) and P* + 1 (odd: x0 = 1/2 on Z + 1/2, o =
-1/2 wrong-signed), S at P_eff = P*. By hand: at large tau Bad runs
to |o| about A_1/(2 pi) (165 at T(P*)) against N about A_1/2,
the window's 4|u| against sin's 2 pi |u|; so the largest ratio in
[0.15, 0.6], past tau = 1; at 1/1229 below 0.01. K_x in mpmath at
60 digits (u near 3e-28).
(C2) PLANT: N := 0 (the window collapsed to o = 0 and the wrong-signed
class) puts the ratio past 1 at some cell.
READ (section_C1, copied from its print): largest truth/S
0.325828 at tau 0.13963 (P*), 0.317892 at 3.5999 (P* + 1), a
plateau near 0.31 from tau 0.14 to 3.6 at both, not a peak past 1
(the place predicted, missed); at the first centre, 0.00082535,
1.095e-4; farthest Bad |o|
153 at both x0 (152.5 at 1/2, P* + 1), none at 300..400; wrong-
signed 2P|K_x|/Gamma 1 at P* + 1 (within 1e-12, a bar set after
this read; the mean value's
first order, which no 60-digit compare splits, so the bound is K7's
F by hand), 4.9e-32 at P* (o = 0's rounding, s phi 0 there); the
plant largest 153.5, past 1 at 85 and 83 of 100 cells. Every bar
and the plant hold: S dominates the truth at P*. The rehearsal
caught the sign's residue mod P dropped at x0 = 0, o < 0, and
Gamma's Fourier series truncated below 0.05 (images now).
(C4) K6's inequalities, mpmath at 40 digits (frozen before code):
(a) |K_x|/(lambda |sin 2 pi x| K) on 40 geometric tau in [1/1229,
20] by 200 x in (0, 1/2) geometric toward both ends: least at least
1, within 1e-6 of it (x near 0, large tau); at tau = 1/1229, x =
1/2 - 1e-6 in [290, 300] (by hand (1/(16 tau^2) - 1/(2 tau))/(2 pi
lambda) = 295.6). PLANT: (1 + q)^3 read as 1 + q puts the least
below 1 at every tau <= 1. AFTER THE FIRST RUN: the least read 1 -
2.3e-41 at tau 9.19: near x = 0 the ratio is 1 + 2q(1 - cos 2 pi x)
+ ..., about 1 + 4e-11 q at x = 1e-6, under 40 digits once q <
2.5e-30 (tau > 1.74), where the least is rounding and its place
moves between runs; so (a) reads at least 1 - 1e-35. Where q shows,
least - 1 is 0.94 at 1/1229 falling to 2.1e-26 at 0.89, each such
row's least at x = 1e-6.
(b) (b/K)/(rho w(1/(4 tau))), images |k| <= 6, K4's Abar as
written, on 30 geometric tau in [1/1229, 1/16] by 400 x in [0,
1/2]: largest at most 1, in [0.30, 0.40]. PLANT: w at 1/(8 tau)
(a = |xi|/(4 tau)) puts it past 1 at every tau.
(c) b(x0 + u) at u = +-u_w against bloc_x0, same tau: at most 1,
in [0.97, 1] at x0 = 1/2, [0.5, 1] at 0. PLANT: the pair's 2
dropped puts x0 = 1/2 past 1.
(d) Printed: N(1/1229) (by hand 0.26 1229^2, in [3.5e5, 4.5e5]),
and the largest u_w/tau over the tau read, at most 4e-22.
(C5) K4's and K5's majorants against W1's and W2's own at P0 = P*
(frozen before code): (a) majorant.image_A (uniform arm, both
signs) over Abar on 60 geometric a in [1e-3, (K_img + 1)c/2] by 60
geometric tau in [1/1229, 1/4]: largest at most 1, in [0.55,
0.572] (|l|/(1.75(1.5a + a^3 tau)) near 1/1.75 where w -> 0 or
w >> 1). PLANT: Abar's a^3 tau dropped puts it past 1. (b) W2's A_k
(majorant.mode_const's constants, mode_A's formula in mpmath at 60
digits) over Abar_k at k = 1..50 (to K_mode(P*) = 48; first run
at 1..40, before K3's modes moved to 1/64) on 60 geometric tau in [1/64,
T(P*)] and at each tau*_k = 3/(8 pi^2 k^2): largest at most 1, at
tau*_k in [0.95, 1]. PLANT: C_k's 2.167(2 pi k)^2 dropped puts
tau*_k's ratio past 1 at every k.
READ (section_C5, copied from its print): (a) a to 14133.5, largest
A/Abar 0.5714285708086326 (doubles read W1's leading order; its
1e-26 terms sit inside 1.75 against 1.7183), the plant 1.902e7;
(b) largest A_k/Abar_k
1.0 (below it in 60 digits), at tau*_k 0.996632 (k = 1) to
0.999664 (k = 50; 0.999646 at k = 40 in the first run), the
plant's least 1.551. Every bar and plant
holds: K4 and K5 stand, K5 within 0.4% of W2 at its corner.
READ (section_C4, copied from its second print): (a) least - 1
-3.4438e-41 at tau 15.433 (the wall), 295.728 at 1/1229 by 1/2 -
1e-6, the plant below 1 at every tau <= 1; (b) largest 0.333333
(b/K is w(1/(4 tau)) at 1/2, rho 3), the plant's least 1.894; (c)
x0 = 0 in [0.55, 1.0], x0 = 1/2 in [0.985891, 1.0], the plant past
1 at every tau; (d) N(1/1229) 384475, u_w/tau 3.529e-22 (the grid's;
the sup is 3.538e-22, K6 (iii)). Every
bar and plant holds: K6 stands.
(C3) The past_T bound (K2), in mpmath at 120 digits, at P = 1e10, P*,
1e40 (frozen before code): (a) against the k = 1..3 terms summed
exactly in mpmath (a lower bound on the truth): bound/sum - 1 in
[0, 1e-9]; (b) against uniform.past_T in doubles: in [1 - 1e-12,
1 + 1e-9], the float wall the read's own rounding, 1.3e-13 past
1e26; (c) eps_T(P*) in [3.5e-28, 4e-28]. PLANT: 3/d dropped from
the k = 1 exponent puts (a) below 0 at all three P.
READ (section_C3, copied from its print): eps_T 2.23657e-8 at 1e10,
3.74746e-28 at P*, 6.38122e-38 at 1e40; bound/sum - 1 2.20657e-18,
2.78233e-58, 6.35122e-78 (the bound second-order tight); bound/read
0.9999999999999845, 1.0000000000000067, 1.0000000000000124 (the
float wall); the plant's bound/sum - 1 -2.20657e-8, -3.72505e-28,
-6.35122e-38. (a), (b), (c) and the plant hold: K2 stands, and P
past_T <= 2.30897e-13 (1 + 3.75e-28) P at every P >= P*.
(C6) K7's C_K (frozen before code; section_CK, mpmath at 30 digits,
2000 geometric cells a band-free span, the cuts added): (a) per
cell, S_sup at least S (S_closed) at t1, the centre and t2 (every
cell); PLANTS, one per arm: e^(-z^2/(4 tau)) read at t1 fails (a)
on some image cell; q^(k^2) read at t2 fails it on some mode
cell. (b) Known answer: the same cells on f0's shape tau^-5.5
e^(-1/(16 tau)) (tau^-5.5 at t1, the exponential at t2) against
16^4.5 (Gamma(4.5, 1/(16 tau1)) - Gamma(4.5, 1/(16 tau0))), the
exact int: ratio in [1, 1.1]. (c) The tails'
majorants at their ends: S(tau0)/f0(tau0) and S(tau1)/f_inf(tau1)
at most 1; PLANT: f0 with 0.557 for 0.557 tau^-2 puts the first
past 1. (d) Printed: C_K (predicted in [1e3, 1e5], the cells near
1/88 most of it), L, both tails, and the cells' sum against the
midpoint sum of S at the centres (at most 1.2 times it).
MARGIN: m(P) = 1/4 - F_K1(P) - 2.06e-4 - (C_K + L/P^2 + 8e-40
T(P) + W(P))/P - L_T (1 + eps_T(P)), F_K1 = (1/400)(1/P + 4/(c(1
- c/P))) e^(c^2/(16 P)) (K1's with its eps kept), on P = P* 10^j
to 1e300: predicted least at P*, in [0.2497, 0.2498].
READ (section_CK, copied from its print): 2002 cells, none
with S_sup below S; the plants fail 535 (e at t1) and 1254 (q at
t2) cells; known answer 1.0285808; S/f0 0.867988 at tau0, S/f_inf
0.95886 at tau1, the plant's S/f0 3.47e6; C_K 7809.09775 (tails
1.4e-44, 2.0e-678), L 1.30067e8, cells/midpoint 1.03648, the
largest cell 68.58 at tau 0.014457, 0.9762 of the cells' sum in
[1/176, 1/44], S's largest at 0.011666; the margin 0.249785863303 at
P*, least there, rising to 0.249793105305 at 1.3e300 (F_K1
8.1367e-6, C_K/P* 5.8e-27, W/P* 7.4e-30). Every bar and plant
holds and the KILL misses: the margin is positive at every P >=
P*.
(C7) The residues (frozen before code). The bound at d = P + 1 needs
floor(H) <= C_d = floor(P/4), and H <= W + E + the larger excess,
so the margin at d is at least b(P) - F - X - R_maj, b = C_d + 1 -
W, W = P/4 at even P and (P^2 - 1)/(4P) + 1/P^2 at odd (the
sawtooth's, classes.py). By hand, at P = 0, 1, 2, 3 mod 4 (d =
1, 2, 3, 0 mod 4): b = 1, 3/4 + 1/(4P) - 1/P^2, 1/2, 1/4 + 1/(4P) - 1/P^2 =
budget(P), the least of the four at every P >= 2. So a margin read
at budget(P), here and in box.py, holds at every d, not only d = 0
mod 4. Read in exact rationals at P = 399..20000 and P* + j, j =
0..3: b equal to its residue's form, at least budget(P), and equal
to it at P = 3 mod 4 alone. PLANT: C_d read as floor((P - 1)/4)
puts b below budget(P) at every P = 0 mod 4.
READ (section_C7, copied from its print): least b - budget by P =
0, 1, 2, 3 mod 4 0.749381, 0.5, 0.249384, 0; the forms and the
bound hold, the plant fires at all 4902 P = 0 mod 4 read.
KILL: the closed form's margin at or below 0 at some P >= P*
(read on a grid of P to 1e300, with each piece's monotonicity in P
proved by hand), or C1's ratio above 1.
FINDINGS. The closed form dominates the true class sum at P* and P* + 1
(largest truth/S 0.326, a plateau near 0.31 from tau 0.14 to 3.6), every
piece K1..K7 holds against its control and its plant fires, and the
margin is positive at every P >= P*, least at P*, 0.2497859. With
box.py below P*, the real support's margin is positive at every d = P +
1 = 0 mod 4 from 400 on, under both scripts' conditions, and so (C7)
at every d from 400 on, the other residues' budgets being larger.
A bare run is the verdict: C3, C4, C5, C1 (with C2), C6, the margin
and C7 in turn, each bar and plant a check; it exits 0 iff every check passes.
RUN RECORD: 36 of 36 checks, every READ above reprinted to its digits;
121.2 s, peak 45.9 MB.
"""
import math
import os
import sys
from fractions import Fraction as Fr
os.environ.setdefault("OPENBLAS_NUM_THREADS", "1")
import mpmath as mp # noqa: E402
import numpy as np # noqa: E402
import majorant as MJ # noqa: E402
import uniform as UN # noqa: E402
from every import T_of # noqa: E402
CHECKS = []
def check(name, ok, detail=""):
CHECKS.append(bool(ok))
print(f" [{'PASS' if ok else 'FAIL'}] {name}" +
(f" -- {detail}" if detail else ""), flush=True)
P_STAR = 399 * 2 ** 25 * 10 ** 20
mp.mp.dps = 120
def L_T():
return 8 * mp.mpf(399) ** 2 * mp.exp(-4 * mp.pi ** 2) / (4 * mp.pi ** 2)
def T_mp(P):
return 1 + mp.log(mp.mpf(P) / 399) / (2 * mp.pi ** 2)
def eps_T(P, plant=False):
"""(K2) past_T(P, T(P)) <= L_T (1 + eps_T(P)); the plant drops 3/d
from the k = 1 exponent."""
P = mp.mpf(P)
d, T, pi2 = P + 1, T_mp(P), mp.pi ** 2
eta = 13 * pi2 / (3 * d * d)
lin = (0 if plant else 3 / d) + eta
k1 = (d / P) ** 3 * mp.exp(4 * pi2 * T * lin) / ((1 - eta)
* (1 - 2 * pi2 / d ** 3))
a = 2 * pi2 * T * (P / d) ** 3
e2 = (((d / P) ** 3 / 2 + 4 * pi2 / P ** 3)
* mp.exp(-4 * pi2 * T * (1 - 6 / d)) / (1 - mp.exp(-4 * a)))
e3 = ((d / 2) * (1 + 2 * d / mp.log(3)) * (4 * pi2 / P ** 3)
* mp.exp(4 * pi2 * T - T * P ** 3 * mp.log(3) / (2 * d)))
lb = min(1.196 * (2 * mp.pi * 2001 / d) ** 2, mp.mpf(0.846))
e4 = ((d / 2) * (1 + 2 * d / lb) * (4 * pi2 / P ** 3)
* mp.exp(4 * pi2 * T - T * P ** 3 * lb / (2 * d)))
return k1 - 1 + e2 + e3 + e4
def past_T_mp(P, kmax=3):
"""The Q-bound's terms k = 1..kmax, exactly, at T = T(P)."""
P = mp.mpf(P)
d, T = P + 1, T_mp(P)
out = mp.mpf(0)
for k in range(1, kmax + 1):
x = mp.log1p(min(8 * mp.sin(mp.pi * k / d) ** 2, 2))
v = x / (2 * d)
out += 8 / P * mp.exp(-T * P ** 3 * x / (2 * d)) / -mp.expm1(-v)
return out
def section_C3():
print("(C3) the past_T bound (K2) against the truth's first terms and "
"uniform.past_T")
rows = []
for P in (10 ** 10, P_STAR, 10 ** 40):
e = eps_T(P)
bound = L_T() * (1 + e)
tr = past_T_mp(P)
rd = UN.past_T(P, T_of(P))
pl = L_T() * (1 + eps_T(P, plant=True))
ra, rb, rp = bound / tr - 1, float(bound) / rd, pl / tr - 1
rows.append((P, e, ra, rb, rp))
print(f" P {float(P):.4e}: eps_T {mp.nstr(e, 6)}, bound/sum - 1 "
f"{mp.nstr(ra, 6)}, bound/read {rb!r}, plant bound/sum - 1 "
f"{mp.nstr(rp, 6)}")
check("C3 (a) bound/sum - 1 in [0, 1e-9] at all three",
all(0 <= r[2] <= 1e-9 for r in rows))
check("C3 (b) bound/read in [1 - 1e-12, 1 + 1e-9] at all three",
all(1 - 1e-12 <= r[3] <= 1 + 1e-9 for r in rows))
check("C3 (c) eps_T(P*) in [3.5e-28, 4e-28]",
3.5e-28 <= rows[1][1] <= 4e-28)
check("C3 PLANT: (a) below 0 at all three", all(r[4] < 0 for r in rows))
def lam(t, plant=False):
"""(K6) lambda(tau); the plant reads (1 + q)^3 as 1 + q."""
q = mp.exp(-4 * mp.pi ** 2 * t)
return 4 * mp.pi * q / ((1 + q) ** (1 if plant else 3) * (1 - q))
def G0(y, t):
return mp.exp(-y * y / (4 * t)) / mp.sqrt(4 * mp.pi * t)
def K_Kx(x, t):
"""K and K_x: Gaussian images |k| <= 6 below tau = 0.05 (the pair at
1/2 cancels in K_x with no loss), the Fourier series above."""
if t < 0.05:
ys = [x + k for k in range(-6, 7)]
return (sum(G0(y, t) for y in ys),
sum(-y / (2 * t) * G0(y, t) for y in ys))
q, K, Kx = mp.exp(-4 * mp.pi ** 2 * t), mp.mpf(1), mp.mpf(0)
for k in range(1, 12):
K += 2 * mp.cos(2 * mp.pi * k * x) * q ** (k * k)
Kx -= 4 * mp.pi * k * mp.sin(2 * mp.pi * k * x) * q ** (k * k)
return K, Kx
def w_img(a, t):
"""(K6 (i)) w(a) = a Abar(a, tau), K4's Abar as written."""
return a * (mp.mpf(1.75) * (1.5 * a + a ** 3 * t) + mp.mpf(10) ** -20)
def rho(t):
return 1 + 2 * sum((2 * j - 1) ** 4 * mp.exp(-j * (j - 1) / (4 * t))
for j in range(1, 21))
def b_img(x, t):
return sum(G0(x + k, t) * w_img(abs(x + k) / (2 * t), t)
for k in range(-6, 7))
def N_img(t):
return rho(t) * w_img(1 / (4 * t), t) / (2 * lam(t))
def u_w(t):
P_eff = max(mp.mpf(P_STAR), mp.exp(1 / (16 * t) - mp.mpf(1) / 8) / 1600)
return (N_img(t) + mp.mpf(1) / 2) / P_eff
def bloc(x0, t, uw, plant=False):
"""(K6 (iii)) sup of b_img within u_w of x0; the plant drops the pair's
2 (one image of it)."""
zs = ([j + mp.mpf(1) / 2 for j in range(20)] if x0 else
[mp.mpf(j) for j in range(1, 21)])
terms = [(2 * z) ** 4 * mp.exp(-z * z / (4 * t)) for z in zs]
if plant and x0:
terms[0] /= 2
s = 2 * w_img(1 / (4 * t), t) * sum(terms)
if not x0:
s += w_img(uw / (2 * t), t)
return (1 + mp.mpf(10) ** -20) * s / mp.sqrt(4 * mp.pi * t)
def section_C4():
print("(C4) K6's inequalities, mpmath at 40 digits")
mp.mp.dps = 40
lo = mp.mpf(1) / 1229
ts = [lo * (20 / lo) ** (mp.mpf(i) / 39) for i in range(40)]
g = [mp.mpf(10) ** (-6 + (mp.log10(0.25) + 6) * i / 99) for i in range(100)]
xs = g + [mp.mpf(1) / 2 - v for v in g]
least, least_pl, at = mp.inf, [], None
for t in ts:
row, row_pl = mp.inf, mp.inf
for x in xs:
K, Kx = K_Kx(x, t)
base = abs(Kx) / (abs(mp.sin(2 * mp.pi * x)) * K)
r, rp = base / lam(t), base / lam(t, plant=True)
if r < least:
least, at = r, (t, x)
row, row_pl = min(row, r), min(row_pl, rp)
least_pl.append((t, row_pl))
K, Kx = K_Kx(mp.mpf(1) / 2 - mp.mpf(10) ** -6, lo)
top = abs(Kx) / (abs(mp.sin(2 * mp.pi * (mp.mpf(1) / 2 - mp.mpf(10) ** -6)))
* K * lam(lo))
print(f" (a) least - 1 {mp.nstr(least - 1, 5)} at tau {mp.nstr(at[0], 5)}"
f", x {mp.nstr(at[1], 5)}; at 1/1229, x = 1/2 - 1e-6: "
f"{mp.nstr(top, 6)}")
check("C4 (a) least in [1 - 1e-35, 1 + 1e-6] (the 40-digit wall)",
1 - mp.mpf(10) ** -35 <= least <= 1 + 1e-6)
check("C4 (a) at 1/1229, x = 1/2 - 1e-6, in [290, 300]",
290 <= top <= 300)
check("C4 (a) PLANT: least below 1 at every tau <= 1",
all(r < 1 for t, r in least_pl if t <= 1))
tb = [lo * ((mp.mpf(1) / 16) / lo) ** (mp.mpf(i) / 29) for i in range(30)]
xb = [mp.mpf(i) / 798 for i in range(400)]
big, pl_rows, c_rows = 0, [], []
for t in tb:
bound = rho(t) * w_img(1 / (4 * t), t)
bound_pl = rho(t) * w_img(1 / (8 * t), t)
rows = [b_img(x, t) / K_Kx(x, t)[0] for x in xb]
big = max(big, max(rows) / bound)
pl_rows.append(max(rows) / bound_pl)
uw = u_w(t)
c_rows.append((t, uw / t,
*(max(b_img(x0 + s * uw, t) for s in (1, -1))
/ bloc(x0, t, uw) for x0 in (mp.mpf(0), mp.mpf(1) / 2)),
b_img(mp.mpf(1) / 2 + uw, t)
/ bloc(mp.mpf(1) / 2, t, uw, plant=True)))
check("C4 (b) largest in [0.30, 0.40]", 0.30 <= big <= 0.40,
mp.nstr(big, 6))
check("C4 (b) PLANT: past 1 at every tau", all(r > 1 for r in pl_rows),
f"least {mp.nstr(min(pl_rows), 4)}")
c0 = [r[2] for r in c_rows]
c12 = [r[3] for r in c_rows]
print(f" (c) x0 = 0 in [{mp.nstr(min(c0), 4)}, {mp.nstr(max(c0), 4)}], "
f"x0 = 1/2 in [{mp.nstr(min(c12), 6)}, {mp.nstr(max(c12), 6)}]")
check("C4 (c) x0 = 1/2 in [0.97, 1]", 0.97 <= min(c12) and max(c12) <= 1)
check("C4 (c) x0 = 0 in [0.5, 1]", 0.5 <= min(c0) and max(c0) <= 1)
check("C4 (c) PLANT: x0 = 1/2 past 1 at every tau",
all(r[4] > 1 for r in c_rows))
n0 = N_img(lo)
uwt = max(r[1] for r in c_rows)
check("C4 (d) N(1/1229) in [3.5e5, 4.5e5]", 3.5e5 <= n0 <= 4.5e5,
mp.nstr(n0, 6))
check("C4 (d) largest u_w/tau at most 4e-22", uwt <= 4e-22,
mp.nstr(uwt, 4))
mp.mp.dps = 120
def K_img(P):
return math.ceil(math.sqrt(4 * math.log(P) + 100)) + 2
def c_of(P):
return math.ceil(16 * math.log(1600 * P)) + 1
def Abar(a, t, plant=False):
"""(K4) the images' majorant at P0 = P; the plant drops a^3 tau."""
return 1.75 * (1.5 * a + (0 if plant else a ** 3 * t)) + 1e-20
def Abar_k(k, t, plant=False):
"""(K5) the modes' majorant at P0 = P, P-free; the plant drops C_k's
2.167 (2 pi k)^2."""
w = 2 * mp.pi * k
b1, b2 = 3 * mp.pi * k, w ** 3 * t
C = (mp.sqrt(b1 * b2 ** 3 / 3) + (0 if plant else 2.167 * w ** 2)
+ 3.06 * w + (2.04 * w ** 3 + 1.084 * w ** 4) * t)
return abs(b1 - b2) + C / P_STAR
def section_C5():
print("(C5) K4's and K5's majorants against W1's and W2's own at P0 = P*")
P0 = float(P_STAR)
amax = (K_img(P0) + 1) * c_of(P0) / 2
a = np.geomspace(1e-3, amax, 60)
big = big_pl = 0.0
for t in np.geomspace(1 / 1229, 0.25, 60):
for s in (1.0, -1.0):
A = MJ.image_A(a, t, s, P0)
big = max(big, float(np.max(A / Abar(a, t))))
big_pl = max(big_pl, float(np.max(A / Abar(a, t, plant=True))))
print(f" (a) a in [1e-3, {amax}]: largest A/Abar {big!r}")
check("C5 (a) largest A/Abar in [0.55, 0.572]", 0.55 <= big <= 0.572)
check("C5 (a) PLANT: past 1", big_pl > 1, f"{big_pl:.4g}")
mp.mp.dps = 60
T0 = mp.mpf(T_of(P0))
grid = [mp.mpf(1) / 64 * (T0 * 64) ** (mp.mpf(i) / 59) for i in range(60)]
Pm = mp.mpf(P_STAR)
worst, star, star_pl = 0, [], []
for k in range(1, 51):
th = 2 * math.pi * k / P0
sg, mb, ec, rho_ = MJ.mode_const(th, P0)
Ec, zc = mp.mpf(P0 * ec), mp.mpf(rho_ * P0 * P0)
b1 = 3 * mp.pi * k
def A_W2(t):
b2 = (2 * mp.pi * k) ** 3 * t
return (mp.sqrt((b1 - b2) ** 2 + b1 * b2 ** 3 / (3 * Pm ** 2))
+ Ec * mp.exp(zc * t) + mp.sqrt(1 + 2.25 * mp.mpf(th) ** 2)
* Pm * mp.expm1(zc * t))
for t in grid:
worst = max(worst, A_W2(t) / Abar_k(k, t))
ts = 3 / (8 * mp.pi ** 2 * k * k)
star.append(A_W2(ts) / Abar_k(k, ts))
star_pl.append(A_W2(ts) / Abar_k(k, ts, plant=True))
worst = max(worst, star[-1])
print(f" (b) k = 1..50: largest A_k/Abar_k {mp.nstr(worst, 8)}; at tau*_k "
f"in [{mp.nstr(min(star), 6)}, {mp.nstr(max(star), 6)}] (k = 1: "
f"{mp.nstr(star[0], 6)}, k = 50: {mp.nstr(star[-1], 6)})")
check("C5 (b) largest A_k/Abar_k at most 1", worst <= 1)
check("C5 (b) at tau*_k in [0.95, 1]",
0.95 <= min(star) and max(star) <= 1)
check("C5 (b) PLANT: past 1 at every k", all(r > 1 for r in star_pl),
f"least {mp.nstr(min(star_pl), 4)}")
mp.mp.dps = 120
def Kxx_half(t):
"""Gamma = |K_xx(1/2, tau)|: images |k| <= 6 below tau = 0.05, the
Fourier series above (as K_Kx)."""
if t < 0.05:
return abs(sum(G0(y, t) * (y * y / (4 * t * t) - 1 / (2 * t))
for y in (mp.mpf(1) / 2 + k for k in range(-6, 7))))
q = mp.exp(-4 * mp.pi ** 2 * t)
return abs(sum(-8 * mp.pi ** 2 * k * k * (-1) ** k * q ** (k * k)
for k in range(1, 12)))
def b_true_img(x0, u, t, P0):
"""(C1) P^3 times the images' arm at x0 + u: W1's A at P0, |k| <=
K_img(P0); the k = 0 image's xi from u, the rest from x0."""
K = K_img(P0)
ks = np.array([k for k in range(-K, K + 1) if x0 + k != 0], float)
xi = np.concatenate([x0 + ks, [u] if x0 == 0 else []])
xi = xi[xi != 0]
a = np.abs(xi) / (2 * t)
G = (4 * np.pi * t) ** -0.5 * np.exp(-xi * xi / (4 * t)) * a
return float((G * np.minimum(MJ.image_A(a, t, np.sign(xi), P0), 1e300)
).sum())
def b_true_mode(t, P0):
"""(C1) P^3 times the modes' arm: W2's A_k at P0, k <= K_mode(P0)."""
if P0 not in MODE_CST:
Km = math.ceil(math.sqrt(7.34 * (3 * math.log(P0) + 95)))
MODE_CST[P0] = []
for k in range(1, Km + 1):
th = 2 * math.pi * k / P0
sg, mb, ec, rho_ = MJ.mode_const(th, P0)
MODE_CST[P0].append((k, th, P0 * ec, rho_ * P0 * P0))
return sum(4 * math.pi * k * math.exp(-4 * math.pi ** 2 * k * k * t)
* float(MJ.mode_A(k, th, Ec, zc, t, P0))
for k, th, Ec, zc in MODE_CST[P0])
MODE_CST = {}
def S_closed(t, plant=False):
"""(K6) S(tau) at P_eff: S_img to 1/16, S_mode from 1/64, the smaller
between; the plant (C2) puts N := 0."""
out = []
if t <= mp.mpf(1) / 16:
N = 0 if plant else N_img(t)
uw = u_w(t)
out.append(2 * (2 * N + 2) * (bloc(0, t, uw) + bloc(1, t, uw)))
if t >= mp.mpf(1) / 64:
q = mp.exp(-4 * mp.pi ** 2 * t)
bm = sum(4 * mp.pi * k * q ** (k * k) * Abar_k(k, t)
for k in range(1, 61))
K_lo = max(2 * G0(mp.mpf(1) / 2, t), 1 - 2 * q)
N = 0 if plant else bm / (2 * lam(t) * K_lo)
out.append(4 * (2 * N + 2) * bm)
return min(out)
def truth(t, P, scan=400):
"""(C1) the true class sum at P: per (x0, o), u = o/P apart from x0,
Bad iff s phi <= B + fl, each Bad class 2(B + fl); with the wrong-
signed classes' 2P|K_x| apart, and the farthest Bad |o| per x0."""
tf = float(t)
fl = 4e-40 / P if tf < 0.25 else 4e-40 * math.exp(-4 * math.pi ** 2
* tf) / P ** 2
bm = b_true_mode(tf, float(P)) if tf >= 1 / 64 else math.inf
tot, wrong, far = 0.0, 0.0, {}
for x0 in (0, 1): # x0 = 0, 1/2
half = x0 and P % 2
os_ = [mp.mpf(j) + (mp.mpf(1) / 2 if half else 0)
for j in range(-scan - (1 if half else 0), scan + 1)]
bi = (b_true_img(0.5, 0.0, tf, float(P)) if x0 and tf < 0.25
else math.inf)
far[x0] = -1.0
for o in os_:
u = o / P
x = mp.mpf(x0) / 2 + u
if x0 == 0 and tf < 0.25:
bi = b_true_img(0.0, float(u), tf, float(P))
B = min(bi, bm)
r = mp.mpf(x0) * P / 2 + o + (P if x0 == 0 and o < 0 else 0)
s = 1 if 2 * r < P - 1 else -1
sp = float(s * (-P * K_Kx(x, t)[1]))
if sp <= B + fl:
tot += 2 * (B + fl)
if sp < 0:
wrong += -2 * sp
far[x0] = max(far[x0], float(abs(o)))
return tot, wrong, far
def section_C1(n_cells=300):
print("(C1) the true class sum at P*, P* + 1 against K6's S; (C2) N := 0")
mp.mp.dps = 60
lo = mp.mpf(1) / c_of(P_STAR)
hi = mp.mpf(T_of(P_STAR))
e = [lo * (hi / lo) ** (mp.mpf(i) / n_cells) for i in range(n_cells + 1)]
ts = [mp.sqrt(e[i] * e[i + 1]) for i in range(0, n_cells, 3)]
rows = []
for P in (P_STAR, P_STAR + 1):
for t in ts:
tr, wr, far = truth(t, P)
S, Sp = S_closed(t), S_closed(t, plant=True)
g = Kxx_half(t)
rows.append((P, float(t), tr / float(S), tr / float(Sp),
float(wr / g) if wr else 0.0, far[0], far[1]))
for P in (P_STAR, P_STAR + 1):
rp = [r for r in rows if r[0] == P]
big = max(rp, key=lambda r: r[2])
print(f" P {'P*' if P == P_STAR else 'P* + 1'}: largest truth/S "
f"{big[2]:.6g} at tau {big[1]:.5g}; at tau {rp[0][1]:.5g} "
f"{rp[0][2]:.4g}; farthest Bad |o| at 0: "
f"{max(r[5] for r in rp):g}, at 1/2: {max(r[6] for r in rp):g}"
f"; wrong-signed 2P|K_x|/Gamma largest "
f"{max(r[4] for r in rp):.6g}")
print(f" plant (N := 0): largest {max(r[3] for r in rp):.4g}, "
f"cells past 1: {sum(r[3] > 1 for r in rp)} of {len(rp)}")
for r in rows[::10]:
print(" ", " ".join(f"{v:.4g}" for v in r))
big = max(r[2] for r in rows)
check("C1 largest truth/S at most 1", big <= 1, f"{big:.6g}")
check("C1 largest truth/S in [0.15, 0.6]", 0.15 <= big <= 0.6)
check("C1 no Bad class at 300 <= |o| <= 400",
all(max(r[5], r[6]) < 300 for r in rows))
check("C1 wrong-signed 2P|K_x|/Gamma within 1e-12 of 1 at P* + 1 (first"
" order; K7 bounds it by F, not Gamma)",
all(abs(r[4] - 1) <= 1e-12 for r in rows if r[0] == P_STAR + 1))
check("C2 PLANT: past 1 at some cell", any(r[3] > 1 for r in rows))
mp.mp.dps = 120
TAU0, TAU1 = mp.mpf(1) / 2000, mp.mpf(40)
def P_eff(t):
return max(mp.mpf(P_STAR), mp.exp(1 / (16 * t) - mp.mpf(1) / 8) / 1600)
def S_sup(t1, t2, plant=None):
"""(K7) sup of S on [t1, t2], each monotone factor at its worse end;
the plants read e^(-z^2/(4 tau)) at t1 ("e") or q^(k^2) at t2
("q")."""
out = []
if t2 <= mp.mpf(1) / 16:
w_hi = w_img(1 / (4 * t1), t1)
N = rho(t2) * w_hi / (2 * lam(t2))
te = t1 if plant == "e" else t2
a_w = (N + mp.mpf(1) / 2) / (2 * t1 * P_eff(t2))
pref = (1 + mp.mpf(10) ** -20) / mp.sqrt(4 * mp.pi * t1)
eh = sum((2 * j + 1) ** 4 * mp.exp(-(j + mp.mpf(1) / 2) ** 2
/ (4 * te)) for j in range(20))
e0 = sum((2 * z) ** 4 * mp.exp(-mp.mpf(z) ** 2 / (4 * te))
for z in range(1, 21))
bh = pref * 2 * w_hi * eh
b0 = pref * (w_img(a_w, t2) + 2 * w_hi * e0)
out.append(2 * (2 * N + 2) * (b0 + bh))
if t1 >= mp.mpf(1) / 64:
tq = t2 if plant == "q" else t1
q1 = mp.exp(-4 * mp.pi ** 2 * tq)
bm = mp.mpf(0)
for k in range(1, 61):
b1 = 3 * mp.pi * k
far = max(abs(b1 - (2 * mp.pi * k) ** 3 * t1),
abs(b1 - (2 * mp.pi * k) ** 3 * t2))
Ck = Abar_k(k, t2) - abs(b1 - (2 * mp.pi * k) ** 3 * t2)
bm += 4 * mp.pi * k * q1 ** (k * k) * (far + Ck)
h = mp.mpf(1) / 2
K_lo = max(min(2 * G0(h, t1), 2 * G0(h, t2)),
1 - 2 * mp.exp(-4 * mp.pi ** 2 * t1))
N = bm / (2 * lam(t2) * K_lo)
out.append(4 * (2 * N + 2) * bm)
return min(out)
def f0(t, plant=False):
"""(K7) the tail's majorant below tau0; the plant drops tau^-2 from
2N + 2's bound."""
return (mp.mpf("0.0046") * (t ** 2 if plant else 1) * t ** -5.5
* mp.exp(-1 / (16 * t)))
def f_inf(t):
return 80500 * t ** 3 * mp.exp(-4 * mp.pi ** 2 * t)
def cells_K(n=2000):
e = [TAU0 * (TAU1 / TAU0) ** (mp.mpf(i) / n) for i in range(n + 1)]
e += [mp.mpf(1) / 64, mp.mpf(1) / 16]
return sorted(set(e))
def c_mp(P):
return mp.ceil(16 * mp.log(1600 * mp.mpf(P))) + 1
def F_K1(P):
P = mp.mpf(P)
c = c_mp(P)
return (1 / P + 4 / (c * (1 - c / P))) * mp.exp(c * c / (16 * P)) / 400
def W_of(P):
P = mp.mpf(P)
return 1 / (2 * mp.sqrt(mp.pi * (1 / c_mp(P) - P ** -2)))
def section_CK():
print("(C6) K7's C_K by cells, the tails by hand; the margin")
mp.mp.dps = 30
e = cells_K()
Sv = {t: S_closed(t) for t in e}
CK = L = mid = known = near = mp.mpf(0)
bad, bad_e, bad_q, top = 0, 0, 0, (0, None)
for t1, t2 in zip(e[:-1], e[1:]):
tc = mp.sqrt(t1 * t2)
Sc = S_closed(tc)
s = S_sup(t1, t2)
w = t2 - t1
CK += w * s
L += s
mid += w * Sc
if w * s > top[0]:
top = (w * s, t1)
bad += s < max(Sv[t1], Sc, Sv[t2])
if t2 <= mp.mpf(1) / 16:
bad_e += S_sup(t1, t2, "e") < max(Sv[t1], Sc, Sv[t2])
if t1 >= mp.mpf(1) / 64:
bad_q += S_sup(t1, t2, "q") < max(Sv[t1], Sc, Sv[t2])
known += w * t1 ** -5.5 * mp.exp(-1 / (16 * t2))
if mp.mpf(1) / 176 <= t1 and t2 <= mp.mpf(1) / 44:
near += w * s
exact = 16 ** mp.mpf(4.5) * (mp.gammainc(4.5, 1 / (16 * TAU1))
- mp.gammainc(4.5, 1 / (16 * TAU0)))
tail0 = mp.mpf("0.0046") * 16 ** mp.mpf(4.5) * mp.gammainc(
4.5, 1 / (16 * TAU0))
a = 4 * mp.pi ** 2
tail1 = 80500 * mp.gammainc(4, a * TAU1) / a ** 4
cells = CK
CK += tail0 + tail1
L += f0(TAU0) + f_inf(TAU1)
r0, r1 = Sv[e[0]] / f0(e[0]), Sv[e[-1]] / f_inf(e[-1])
r0p = Sv[e[0]] / f0(e[0], plant=True)
print(f" {len(e) - 1} cells on [{mp.nstr(TAU0, 4)}, {mp.nstr(TAU1, 4)}]"
f"; (a) cells with S_sup below S: {bad}; PLANTS: e at t1 {bad_e}"
f", q at t2 {bad_q}")
print(f" (b) known answer: cells/exact {mp.nstr(known / exact, 8)}")
print(f" (c) S/f0 at tau0 {mp.nstr(r0, 6)}, S/f_inf at tau1 "
f"{mp.nstr(r1, 6)}; PLANT S/f0 {mp.nstr(r0p, 6)}")
print(f" (d) C_K {mp.nstr(CK, 10)} (cells {mp.nstr(cells, 10)}, tails "
f"{mp.nstr(tail0, 4)}, {mp.nstr(tail1, 4)}), L {mp.nstr(L, 6)}, "
f"cells/midpoint {mp.nstr(cells / mid, 6)}, largest cell "
f"{mp.nstr(top[0], 4)} at tau {mp.nstr(top[1], 5)}; share in "
f"[1/176, 1/44] {mp.nstr(near / cells, 4)}, S's largest at tau "
f"{mp.nstr(max(e, key=lambda t: Sv[t]), 5)}")
check("C6 (a) every cell's S_sup at least S", bad == 0)
check("C6 (a) PLANTS fail some cell, each arm", bad_e > 0 and bad_q > 0)
check("C6 (b) known answer in [1, 1.1]", 1 <= known / exact <= 1.1)
check("C6 (c) both tails' ratios at most 1", r0 <= 1 and r1 <= 1)
check("C6 (c) PLANT: past 1", r0p > 1)
check("C6 (d) C_K in [1e3, 1e5]", 1e3 <= CK <= 1e5)
check("C6 (d) cells/midpoint at most 1.2", cells / mid <= 1.2)
mp.mp.dps = 120
eT = eps_T(P_STAR) # eps_T falls in P (K2)
rows = []
for j in range(0, 271):
P = mp.mpf(P_STAR) * mp.mpf(10) ** j
m = (mp.mpf(1) / 4 - F_K1(P) - mp.mpf("2.06e-4")
- (CK + L / P ** 2 + mp.mpf("8e-40") * T_mp(P) + W_of(P)) / P
- mp.mpf("2.30897e-13") * (1 + eT))
rows.append((j, m))
least = min(rows, key=lambda r: r[1])
mono = all(b[1] >= a[1] for a, b in zip(rows, rows[1:]))
print(f" MARGIN: at P* {mp.nstr(rows[0][1], 12)}, least "
f"{mp.nstr(least[1], 12)} at P* 10^{least[0]}, at 1.3e300 "
f"{mp.nstr(rows[-1][1], 12)}; rising along the grid: {mono}")
print(f" pieces at P*: F_K1 {mp.nstr(F_K1(P_STAR), 6)}, C_K/P* "
f"{mp.nstr(CK / P_STAR, 4)}, W/P* {mp.nstr(W_of(P_STAR) / P_STAR, 4)}"
f", eps_T {mp.nstr(eT, 4)}")
check("MARGIN least at P*, in [0.2497, 0.2498]", least[0] == 0
and mp.mpf("0.2497") <= least[1] <= mp.mpf("0.2498"))
check("KILL not fired: the margin above 0 at every P read", least[1] > 0)
def section_C7():
print("(C7) the residues: b = C_d + 1 - W against budget(P)")
def W(P):
return Fr(P, 4) if P % 2 == 0 else Fr(P * P - 1, 4 * P) + Fr(1, P * P)
def bud(P):
return Fr(1, 4) + Fr(1, 4 * P) - Fr(1, P * P)
form = {0: lambda P: Fr(1), 2: lambda P: Fr(1, 2), 3: bud,
1: lambda P: Fr(3, 4) + Fr(1, 4 * P) - Fr(1, P * P)}
Ps = list(range(399, 20001)) + [P_STAR + j for j in range(4)]
ok_form = ok_ge = ok_eq = True
least = {r: None for r in range(4)}
plant = [P for P in Ps if P % 4 == 0
and (P - 1) // 4 + 1 - W(P) < bud(P)]
for P in Ps:
b = P // 4 + 1 - W(P)
ok_form &= b == form[P % 4](P)
ok_ge &= b >= bud(P)
ok_eq &= (b == bud(P)) == (P % 4 == 3)
r = P % 4
if least[r] is None or b - bud(P) < least[r]:
least[r] = b - bud(P)
print(" least b - budget by P mod 4 = 0, 1, 2, 3: " + ", ".join(
f"{float(least[r]):.6g}" for r in range(4)))
check("C7 b equal to its residue's form", ok_form)
check("C7 b at least budget(P), equal at P = 3 mod 4 alone",
ok_ge and ok_eq)
check("C7 PLANT: b below budget at every P = 0 mod 4",
len(plant) == sum(P % 4 == 0 for P in Ps), f"{len(plant)} P")
if __name__ == "__main__":
section_C3()
section_C4()
section_C5()
section_C1()
section_CK()
section_C7()
print(f"{sum(CHECKS)}/{len(CHECKS)} checks passed")
sys.exit(0 if all(CHECKS) else 1)