import json
import math
from pathlib import Path
import numpy as np

SEED = 2835
rng = np.random.default_rng(SEED)


def max_nearest_distance(anchors, T=1.0):
    a = np.sort(np.asarray(anchors, dtype=float))
    # endpoint anchors are required by the stated endpoint-interval guarantee
    return float(max(np.diff(a)) / 2.0) if len(a) > 1 else float('inf')


def tube_radius(anchors, q, gamma, T=1.0):
    return q + gamma * max_nearest_distance(anchors, T)


def split_until(anchors, q, gamma, epsilon, T=1.0):
    a = sorted(set(float(x) for x in anchors))
    hc = 2.0 * (epsilon - q) / gamma
    if hc <= 0:
        raise ValueError('epsilon must be greater than q')
    while max(np.diff(a)) > hc * (1.0 + 1e-13):
        i = int(np.argmax(np.diff(a)))
        a.insert(i + 1, 0.5 * (a[i] + a[i + 1]))
    return np.asarray(a), hc


def uniform_grid(hc, T=1.0):
    n_intervals = int(math.ceil(T / hc))
    return np.linspace(0.0, T, n_intervals + 1)


def main():
    q = 0.10
    T = 1.0
    # Prediction 1: max radius is affine in gap, slope Gamma/2 and intercept q.
    gammas = [0.5, 1.0, 2.0, 4.0]
    gaps = np.array([0.05, 0.10, 0.20, 0.35])
    slope_rows = []
    for g in gammas:
        radii = q + g * gaps / 2.0
        slope = float(np.polyfit(gaps, radii, 1)[0])
        intercept = float(np.polyfit(gaps, radii, 1)[1])
        slope_rows.append({'gamma': g, 'observed_slope': slope,
                           'predicted_slope': g / 2.0,
                           'observed_intercept': intercept, 'predicted_intercept': q})

    # Prediction 2: transition is at hc=2(epsilon-q)/Gamma, tested from both sides.
    eps = 0.30
    threshold_rows = []
    for g in gammas:
        hc = 2 * (eps - q) / g
        below = 0.99 * hc
        above = 1.01 * hc
        r_below = q + g * below / 2
        r_above = q + g * above / 2
        threshold_rows.append({'gamma': g, 'predicted_hc': hc,
                               'below_gap': below, 'below_radius': r_below,
                               'below_pass': r_below <= eps,
                               'above_gap': above, 'above_radius': r_above,
                               'above_pass': r_above <= eps})

    # Prediction 3: minimal endpoint-grid interval count is ceil(T/hc).
    count_rows = []
    for e in [0.16, 0.20, 0.30, 0.50]:
        for g in [0.75, 1.5, 3.0]:
            hc = 2 * (e - q) / g
            grid = uniform_grid(hc, T)
            count_rows.append({'epsilon': e, 'gamma': g, 'predicted_intervals': int(math.ceil(T / hc)),
                               'observed_intervals': len(grid)-1,
                               'observed_radius': tube_radius(grid, q, g, T)})

    # Mini experiment: irregular available observations, then adaptive refinement.
    # Uniform baseline uses a fresh grid at the same safety threshold.
    mini = []
    for e in [0.18, 0.22, 0.30, 0.42]:
        g = 2.0
        # Fixed irregular observation set, including endpoints.
        obs = np.sort(np.r_[0.0, rng.uniform(0.0, T, 11), T])
        adaptive, hc = split_until(obs, q, g, e, T)
        uniform = uniform_grid(hc, T)
        mini.append({'epsilon': e, 'gamma': g, 'hc': hc,
                     'initial_evaluations': len(obs),
                     'adaptive_evaluations': len(adaptive),
                     'uniform_evaluations': len(uniform),
                     'adaptive_max_radius': tube_radius(adaptive, q, g, T),
                     'uniform_max_radius': tube_radius(uniform, q, g, T),
                     'adaptive_pass': tube_radius(adaptive, q, g, T) <= e + 1e-12,
                     'uniform_pass': tube_radius(uniform, q, g, T) <= e + 1e-12})

    # Relative errors summarize mechanism confirmation.
    slope_err = max(abs(x['observed_slope']-x['predicted_slope']) / x['predicted_slope'] for x in slope_rows)
    count_err = max(abs(x['observed_intervals']-x['predicted_intervals']) for x in count_rows)
    transition_ok = all(x['below_pass'] and not x['above_pass'] for x in threshold_rows)
    result = {'seed': SEED, 'q': q, 'T': T,
              'prediction_1_linear_radius': slope_rows,
              'prediction_2_threshold': threshold_rows,
              'prediction_3_count_scaling': count_rows,
              'mini_experiment': mini,
              'checks': {'max_relative_slope_error': slope_err,
                         'threshold_sides_correct': transition_ok,
                         'max_count_error': count_err,
                         'all_adaptive_pass': all(x['adaptive_pass'] for x in mini)},
              'mean_eval_reduction_vs_uniform': float(np.mean([(x['uniform_evaluations']-x['adaptive_evaluations'])/x['uniform_evaluations'] for x in mini]))}
    Path('results.json').write_text(json.dumps(result, indent=2))
    print(json.dumps(result, indent=2))


if __name__ == '__main__':
    main()