import json import math from pathlib import Path import numpy as np # Adaptive reset neural-ODE MVP. The scalar flow is deliberately chosen so that # its exact solution makes the paper's stopping rule quantitatively testable. SEED = 2093 rng = np.random.default_rng(SEED) def q95(x): return float(np.quantile(np.asarray(x), 0.95)) def constant_flow_error(s, lam, lam_hat, y0): return np.abs(np.exp(lam * s) - np.exp(lam_hat * s)) * np.abs(y0) def analytic_crossing(eps, lam, lam_hat, yq, horizon=10.0): # First crossing of the exact supremum/quantile error. Bisection is used # only for the analytic verification; the implementation uses a grid. def f(s): return abs(math.exp(lam*s) - math.exp(lam_hat*s)) * yq - eps if f(horizon) <= 0: return horizon lo, hi = 0.0, horizon for _ in range(80): mid = (lo + hi) / 2 if f(mid) > 0: hi = mid else: lo = mid return (lo + hi) / 2 def grid_crossing(eps, lam, lam_hat, y0, dt=0.002, horizon=10.0): ts = np.arange(dt, horizon + 0.5 * dt, dt) e = np.array([q95(constant_flow_error(t, lam, lam_hat, y0)) for t in ts]) hit = np.flatnonzero(e > eps) return float(ts[hit[0]]) if len(hit) else horizon def verify_math(): y0 = np.linspace(0.8, 1.2, 101) yq = np.quantile(y0, .95) rows = [] # Prediction 1: grid boundary converges to the analytic first crossing. for eps in (.01, .03, .08): exact = analytic_crossing(eps, .5, .35, yq) observed = grid_crossing(eps, .5, .35, y0) rows.append({"kind": "tolerance", "epsilon": eps, "predicted": exact, "observed": observed, "abs_error": abs(exact-observed)}) # Prediction 2: increasing mismatch causes an earlier boundary. mismatch_rows = [] for mh in (.49, .45, .40, .30): exact = analytic_crossing(.03, .5, mh, yq) observed = grid_crossing(.03, .5, mh, y0) mismatch_rows.append({"lambda_hat": mh, "predicted": exact, "observed": observed}) # Prediction 3: at fixed mismatch, tolerance monotonically increases window. tol_exact = [analytic_crossing(e, .5, .35, yq) for e in (.005,.01,.02,.04,.08)] tol_obs = [grid_crossing(e, .5, .35, y0) for e in (.005,.01,.02,.04,.08)] # Strict checks are against the discrete-grid resolution, not an arbitrary win. max_grid_error = max(r["abs_error"] for r in rows) monotone_mismatch = all(mismatch_rows[i]["observed"] > mismatch_rows[i+1]["observed"] for i in range(len(mismatch_rows)-1)) monotone_tol = all(tol_obs[i] < tol_obs[i+1] for i in range(len(tol_obs)-1)) passed = max_grid_error <= .00201 and monotone_mismatch and monotone_tol return {"passed": passed, "y95": float(yq), "tolerance_sweep": rows, "mismatch_sweep": mismatch_rows, "tolerance_monotonicity": {"predicted": tol_exact, "observed": tol_obs}, "max_grid_abs_error": max_grid_error} def teacher_a(t, difficulty=1.0): # A slowly changing vector-field coefficient: the local model is theta*y. return .35 + difficulty * (.30*np.sin(.9*t) + .22*np.sin(2.1*t + .4)) def teacher_solution(times, y0, difficulty=1.0): a = teacher_a(times, difficulty) integ = np.zeros_like(times) integ[1:] = np.cumsum(.5*(a[1:]+a[:-1])*np.diff(times)) return y0[:, None] * np.exp(integ[None, :]) def fit_local_theta(times, ys, difficulty=1.0): # Teacher-forced derivative regression for f_theta(y)=theta*y. a = teacher_a(times, difficulty) y = ys[:, :-1] dy = a[:-1][None, :] * y return float(np.sum(y*dy) / (np.sum(y*y) + 1e-12)) def adaptive_windows(times, truth, y0, eps=.035, min_steps=4, cap_steps=35, difficulty=1.0): windows, thetas = [], [] start = 0 n = len(times)-1 while start < n: stop_cap = min(n, start + cap_steps) # Warm-start analogue: use teacher derivative data in the candidate cap. theta = fit_local_theta(times[start:stop_cap+1], truth[:, start:stop_cap+1], difficulty) offsets = np.arange(min_steps, stop_cap-start+1) chosen = stop_cap for off in offsets: s = times[start+off] - times[start] pred = truth[:, start:start+1] * np.exp(theta*s) err = q95(np.abs(truth[:, start+off] - pred[:, 0])) if err > eps: chosen = start + int(off) break # enforce min window even if the first grid point crosses chosen = max(chosen, min(n, start + min_steps)) windows.append((start, chosen)) thetas.append(theta) start = chosen return windows, thetas def rollout(times, y0, windows, thetas): out = np.empty((len(y0), len(times))) out[:, 0] = y0 for (start, stop), theta in zip(windows, thetas): # At deployment the reset is the model's terminal state, not teacher data. for j in range(start+1, stop+1): out[:, j] = out[:, start] * np.exp(theta*(times[j]-times[start])) return out def mini_experiment(): y0 = np.linspace(.8, 1.2, 128) times = np.linspace(0, 12, 1201) truth = teacher_solution(times, y0) # Shared field is the standard one-field teacher-forced derivative fit. theta = fit_local_theta(times, truth) shared = y0[:, None] * np.exp(theta*times[None, :]) windows, thetas = adaptive_windows(times, truth, y0) adaptive = rollout(times, y0, windows, thetas) shared_rmse = float(np.sqrt(np.mean((shared-truth)**2))) adaptive_rmse = float(np.sqrt(np.mean((adaptive-truth)**2))) shared_final = float(np.sqrt(np.mean((shared[:,-1]-truth[:,-1])**2))) adaptive_final = float(np.sqrt(np.mean((adaptive[:,-1]-truth[:,-1])**2))) boundary_errors = [q95(np.abs(truth[:,s]- (truth[:,s-1] if False else adaptive[:,s]))) for s,_ in windows[1:]] return {"shared_theta": theta, "windows": len(windows), "window_lengths": [b-a for a,b in windows], "shared_rmse": shared_rmse, "adaptive_rmse": adaptive_rmse, "shared_final_rmse": shared_final, "adaptive_final_rmse": adaptive_final, "adaptive_boundary_count": len(boundary_errors), "adaptive_boundary_terminal_error_q95": boundary_errors[:8]} def main(): result = {"seed": SEED, "math_verification": verify_math(), "mini_experiment": mini_experiment()} Path("results.json").write_text(json.dumps(result, indent=2)) print(json.dumps(result, indent=2)) if __name__ == "__main__": main()