1import json
  2import math
  3from pathlib import Path
  4import numpy as np
  5
  6SEED = 2734
  7rng = np.random.default_rng(SEED)
  8
  9# 1-D specialization of a 2-D obstacle margin: Phi(x)=distance to obstacle.
 10# The plant moves toward the obstacle with constant speed v, so dPhi/dt=-v exactly.
 11def certificate_sweep():
 12    gamma_true = 1.0
 13    gamma_est = 1.25  # conservative bound
 14    delta = 0.05
 15    distances = np.linspace(0.15, 2.0, 80)
 16    latencies = np.linspace(0.0, 2.0, 81)
 17    permitted = []
 18    unsafe_permitted = []
 19    for d in distances:
 20        for L in latencies:
 21            contract_permits = L <= d / gamma_est - delta
 22            physically_safe = L <= d / gamma_true + 1e-12
 23            permitted.append(contract_permits)
 24            unsafe_permitted.append(contract_permits and not physically_safe)
 25    # The implication is checked over a dense grid, including the boundary.
 26    conservative_violations = int(np.sum(unsafe_permitted))
 27
 28    # With an intentionally underestimated gamma, the observed permit/fail
 29    # transition is compared to the exact physical boundary d/gamma_true-delta.
 30    under_gamma = 0.75
 31    transitions = []
 32    for d in distances:
 33        ls = np.linspace(0, 2.5, 10001)
 34        permits = ls <= d / under_gamma - delta
 35        unsafe = ls > d / gamma_true
 36        bad = np.where(permits & unsafe)[0]
 37        transitions.append(float(ls[bad[0]]) if len(bad) else float("nan"))
 38    valid = np.isfinite(transitions)
 39    observed_underestimate_boundary = float(np.nanmedian(np.array(transitions)[valid]))
 40    # For this sweep the onset is expected near the smallest boundary where
 41    # under-estimation starts to allow an unsafe handoff: d/gamma_true.
 42    expected_underestimate_boundary = float(np.median(distances[distances / under_gamma - delta > distances / gamma_true]))
 43    return {
 44        "grid_points": len(distances) * len(latencies),
 45        "conservative_contract_violations": conservative_violations,
 46        "conservative_violation_rate": conservative_violations / len(permitted),
 47        "underestimated_gamma": under_gamma,
 48        "underestimate_unsafe_onset_median_L": observed_underestimate_boundary,
 49        "expected_physical_boundary_median_L": expected_underestimate_boundary,
 50    }
 51
 52def derivative_verification():
 53    # Phi(t)=phi0-v*t and gamma_true=v: finite differences satisfy
 54    # dPhi/dt=-gamma_true, hence Phi(t)-Phi(0)+gamma_true*t == 0.
 55    phi0, gamma_true = 1.7, 1.3
 56    ts = np.linspace(0.0, 1.0, 1001)
 57    phi = phi0 - gamma_true * ts
 58    finite_slopes = np.diff(phi) / np.diff(ts)
 59    residual = phi - (phi0 - gamma_true * ts)
 60    return {"max_abs_slope_error": float(np.max(np.abs(finite_slopes + gamma_true))),
 61            "max_contract_residual": float(np.max(np.abs(residual))),
 62            "predicted_unit_rate": 1.0,
 63            "observed_normalized_slope": float(np.max(-finite_slopes / gamma_true))}
 64
 65def gamma_boundary_sweep():
 66    # Prediction: for a fixed latency L and margin delta, permission changes at
 67    # gamma_est = d/(L+delta). Conservative gamma has no unsafe permissions;
 68    # underestimated gamma creates unsafe permissions.
 69    d = 1.0
 70    L = 1.05
 71    delta = 0.05
 72    gamma_true = 1.0
 73    gammas = np.linspace(0.4, 1.8, 141)
 74    unsafe = []
 75    permitted = []
 76    for g in gammas:
 77        p = L <= d / g - delta
 78        permitted.append(p)
 79        unsafe.append(p and L > d / gamma_true)
 80    transition = gammas[np.where(np.array(permitted) == False)[0][0]]
 81    predicted = d / (L + delta)
 82    return {
 83        "predicted_gamma_transition": predicted,
 84        "observed_gamma_transition_grid": float(transition),
 85        "transition_abs_error": abs(float(transition) - predicted),
 86        "unsafe_permissions_at_conservative_gamma": int(sum(unsafe[i] for i,g in enumerate(gammas) if g >= gamma_true)),
 87        "unsafe_permissions_at_underestimated_gamma": int(sum(unsafe[i] for i,g in enumerate(gammas) if g < gamma_true)),
 88        "predicted_unsafe_onset_gamma": predicted,
 89        "observed_first_unsafe_gamma": float(gammas[np.where(np.array(unsafe))[0][-1]]) if np.any(unsafe) else None,
 90    }
 91
 92def evaluation_scaling():
 93    # Prediction: a constant-distance contract with usable interval
 94    # B=d/gamma-delta-L has evaluation rate 1/B, and rate ratios scale as
 95    # the inverse certified budget.
 96    d, delta, L, horizon = 1.0, 0.05, 0.05, 100.0
 97    gammas = np.array([0.5, 1.0, 1.5, 2.0])
 98    rows = []
 99    for g in gammas:
100        B = d / g - delta - L
101        observed = math.ceil(horizon / B)
102        predicted = horizon / B
103        rows.append({"gamma": float(g), "budget": B, "observed_evals": observed,
104                     "predicted_evals": predicted, "rate_error": abs(observed-predicted)/predicted})
105    # Relative scaling is the important quantity and is not affected by ceil.
106    return {"rows": rows, "rate_ratio_observed_gamma2_over_gamma1": rows[3]["observed_evals"] / rows[1]["observed_evals"],
107            "rate_ratio_predicted_gamma2_over_gamma1": rows[3]["predicted_evals"] / rows[1]["predicted_evals"]}
108
109def scheduler_comparison():
110    # Same finite approach episodes, random initial distance and speed. A
111    # handoff is unsafe if latency exceeds the distance remaining at launch.
112    n = 3000
113    d0 = rng.uniform(0.5, 2.0, n)
114    v = rng.uniform(0.5, 1.5, n)
115    L = rng.uniform(0.0, 0.8, n)
116    delta = 0.05
117    gamma = 1.6  # safely above the sampled vmax
118    # ordinary periodic inference: fixed interval, no certificate fallback
119    periodic_h = 0.50
120    periodic_unsafe = L > np.maximum(d0 - v * periodic_h, 0) / v
121    # contract: only launch when the certified remaining interval covers latency;
122    # otherwise fallback/urgent refresh is used, so unsafe handoffs are rejected.
123    budget = d0 / gamma - delta
124    contract_launch = L <= budget
125    contract_unsafe = contract_launch & (L > d0 / v)
126    # Naive event trigger uses a fixed distance threshold, not a rate certificate.
127    naive_launch = d0 > 0.25
128    naive_unsafe = naive_launch & (L > d0 / v)
129    return {
130        "episodes": n,
131        "periodic_unsafe_handoff_rate": float(periodic_unsafe.mean()),
132        "naive_event_unsafe_handoff_rate": float(naive_unsafe.mean()),
133        "contract_unsafe_handoff_rate": float(contract_unsafe.mean()),
134        "contract_launch_fraction": float(contract_launch.mean()),
135        "contract_fallback_fraction": float((~contract_launch).mean()),
136        "periodic_evaluations_per_episode": 1.0,
137        "contract_evaluations_per_episode": float(contract_launch.mean()),
138    }
139
140def main():
141    out = {"seed": SEED, "derivative_verification": derivative_verification(),
142           "math_verification": certificate_sweep(),
143           "gamma_boundary_sweep": gamma_boundary_sweep(),
144           "evaluation_scaling": evaluation_scaling(),
145           "scheduler_comparison": scheduler_comparison()}
146    Path("results.json").write_text(json.dumps(out, indent=2))
147    print(json.dumps(out, indent=2))
148
149if __name__ == "__main__":
150    main()