"""Check T618's index-one wording repair without changing any numbers.

Download generate.py, solver.py, repair-input.json and repair-skill-evidence.json
alongside this file, then run: sage -python repair_checks.py
The input is a before/after document snapshot; no key or live write is needed.
"""
import ast
import json
from decimal import Decimal
from pathlib import Path
import numberdb.sage as numberdb
from sage.symbolic.ring import SR
import sage.functions.all
import sage.calculus.calculus
from sage.functions.trig import sin
from sage.symbolic.constants import pi
from sage.rings.rational_field import QQ
from sage.rings.real_arb import RealBallField
from sage.version import version as sage_version
from generate import LaneEmdenMassRadius

work = Path(__file__).resolve().parent
data = json.loads((work / 'repair-input.json').read_text())
before, after = data['before'], data['after']
evidence = json.loads((work / 'repair-skill-evidence.json').read_text())
expected = before['Numbers']['1']['comment'].replace('linear profile', 'index-one profile')
assert after['Numbers']['1']['comment'] == expected
restored = json.loads(json.dumps(after))
restored['Numbers']['1']['comment'] = before['Numbers']['1']['comment']
assert restored == before
old_source, new_source = data['old_generator'], (work / 'generate.py').read_text()
assert new_source.replace(expected, before['Numbers']['1']['comment']) == old_source
assert ast.dump(ast.parse(old_source)) == ast.dump(ast.parse(new_source.replace(expected, before['Numbers']['1']['comment'])))

xi = SR.var('xi')
theta = sin(xi) / xi
residual = (xi**2 * theta.diff(xi)).diff(xi) + xi**2 * theta
assert residual.simplify_full() == 0
assert theta.limit(xi=0) == 1
assert theta.diff(xi).limit(xi=0) == 0
# Nonzero second derivative disproves linearity in the radial variable.
second = theta.diff(xi, 2).subs(xi=pi).simplify_full()
assert (second - 2/pi**2).simplify_full() == 0
mu = (-xi**2 * theta.diff(xi)).subs(xi=pi).simplify_full()
assert mu == pi
coefficient = (4*pi) / 2 * pi**(-2)
assert (coefficient - 2/pi).simplify_full() == 0

B = RealBallField(500)
computed = B(2) / B.pi()
stored = Decimal(before['Numbers']['1']['number'])
num, den = stored.as_integer_ratio()
midpoint = QQ(num) / QQ(den)
ulp = QQ(10) ** int(stored.as_tuple().exponent)
assert computed.is_finite() and B(computed.rad()) < B(10)**(-140)
assert B(midpoint - ulp) < computed < B(midpoint + ulp)

generator = LaneEmdenMassRadius()
enumeration = list(generator.enumerate())
assert enumeration == [{'n': n} for n in before['Numbers']]
generated = generator.value({'n': '1'}, 100)
assert generated['comment'] == expected
value = generated['number']
assert value.is_finite() and B(value.rad()) < B(10)**(-120)
assert B(midpoint - ulp) < value < B(midpoint + ulp)
assert value.overlaps(computed)

print(json.dumps({
    'table': 'T618', 'checks_passed': True, 'sage_version': sage_version,
    'skill_manifest': evidence['manifest'], 'skill_versions': evidence['versions'],
    'scope': 'Index-one analytic claim and comment only; no other ODE integrations rerun.',
    'profile': 'sin(xi)/xi', 'lane_emden_residual': '0',
    'central_value': '1', 'central_derivative': '0',
    'second_derivative_at_pi': '2/pi^2', 'mass_at_pi': 'pi',
    'coefficient': '2/pi', 'independent_ball_bits': 500,
    'independent_ball': str(computed), 'stored_ulp': str(ulp),
    'difference_in_stored_ulps': str((computed - B(midpoint))/B(ulp)),
    'stored_interval_contains_independent_ball': True,
    'stored_interval_contains_generator_ball': True,
    'ordered_enumeration': enumeration, 'numbers_and_parameters_unchanged': True,
    'generator_changed_only_returned_comment': True,
    'revised_comment': expected,
    'source': 'https://www.astro.princeton.edu/~gk/A403/polytrop.pdf (poly.6), (poly.9), (poly.11b)'
}, indent=2))
