"""Small analytic tangent-branch evasion layer and reproducible toy checks.""" from __future__ import annotations import json from dataclasses import dataclass from pathlib import Path import numpy as np EPS = 1e-12 def rot(theta: float) -> np.ndarray: c, s = np.cos(theta), np.sin(theta) return np.array([[c, -s], [s, c]], dtype=float) def tangent_point(center: np.ndarray, target: np.ndarray, radius: float, side: int): """Return tangent point on the target-facing tangent construction. side=+1/-1 selects the two orientations. Invalid disks (d<=r) return None. """ center, target = np.asarray(center, float), np.asarray(target, float) z = center - target d = float(np.linalg.norm(z)) if d <= radius + EPS: return None alpha = np.arcsin(np.clip(radius / d, -1.0, 1.0)) # Formula in the prompt: R(alpha)^(-q), q=side. # The idea specifies R(alpha)^(-q), giving the two tangent orientations. S = target + (np.sqrt(d*d-radius*radius)/d) * rot(-side * alpha) @ z return S def signed_angle(a: np.ndarray, b: np.ndarray) -> float: return float(np.arctan2(np.cross(a, b), np.dot(a, b))) @dataclass class Candidate: side: int point: np.ndarray turn_time: float straight_time: float time: float def plan(position, heading, target, center, radius, speed, turn_rate): """Create both branches and select the minimum estimated completion time.""" p, h, target, center = map(lambda x: np.asarray(x, float), (position, heading, target, center)) out = [] for side in (-1, 1): S = tangent_point(center, target, radius, side) if S is None: continue turn = abs(signed_angle(h, S-p)) / max(turn_rate, EPS) straight = np.linalg.norm(target-S) / max(speed, EPS) out.append(Candidate(side, S, turn, straight, turn+straight)) selected = min(out, key=lambda x: x.time) if out else None return out, selected def softmin(times, tau): t = np.asarray(times, float) m = np.min(t) return float(m - tau*np.log(np.sum(np.exp(-(t-m)/tau)))) def project_velocity(position, velocity, centers, radius, margin=0.0, lookahead_dt=0.0): """Project inward components; lookahead_dt avoids discrete-step crossing.""" x, v = np.asarray(position, float), np.asarray(velocity, float).copy() for c in np.asarray(centers, float): nvec = x-c dist = np.linalg.norm(nvec) if dist < EPS: continue n = nvec/dist # outward normal # A simulator takes a finite step. Activate slightly outside the disk # if the current velocity could cross it in one step. safe_dist = radius + margin + lookahead_dt*np.linalg.norm(v) if dist <= safe_dist and np.dot(v, n) < 0: v -= np.dot(v, n)*n return v def integrate(position, velocity_fn, centers, radius, dt=0.005, steps=800): x = np.asarray(position, float).copy() min_clearance = np.inf for _ in range(steps): v = velocity_fn(x) x += dt*v min_clearance = min(min_clearance, min(np.linalg.norm(x-c) for c in centers)-radius) return x, float(min_clearance) def main(): rng = np.random.default_rng(1218) # Prediction 1: tangent target distance is sqrt(d^2-r^2). rel_errors, clear_errors = [], [] for _ in range(1000): d = rng.uniform(1.01, 10.0); r = rng.uniform(.05, .8) if r >= d: r = .5*d ang = rng.uniform(-np.pi, np.pi) target = rng.normal(size=2); center = target + d*np.array([np.cos(ang), np.sin(ang)]) expected = np.sqrt(d*d-r*r) for side in (-1,1): S = tangent_point(center,target,r,side) rel_errors.append(abs(np.linalg.norm(S-target)-expected)/expected) clear_errors.append(abs(np.linalg.norm(S-center)-r)) # Prediction 2: branch equality occurs when heading points along the # angular bisector of the two tangent rays. For a symmetric disk this is # heading angle zero; sweep headings and locate the sign change. pos=np.array([0.,0.]); target=np.array([10.,0.]); center=np.array([5.,0.]) headings=np.linspace(-.8,.8,1601) diffs=[] for a in headings: cs,_=plan(pos,np.array([np.cos(a),np.sin(a)]),target,center,1.,1.,1.) # side - minus side +; at heading zero both are equal diffs.append(cs[0].time-cs[1].time) diffs=np.asarray(diffs) switch=headings[np.argmin(abs(diffs))] # Prediction 3: projection prevents boundary penetration independent of # inward residual size, modulo Euler integration error. c=np.array([[0.,0.]]) raw_clear=[]; projected_clear=[] x0=np.array([1.0,0.0]); r=.8 for mag in np.linspace(.1,8.,12): raw=integrate(x0,lambda x,m=mag: np.array([-m,0.]),c,r) proj=integrate(x0,lambda x,m=mag: project_velocity(x,np.array([-m,0.]),c,r,lookahead_dt=.005),c,r) raw_clear.append(raw[1]); projected_clear.append(proj[1]) # Prediction 3 (paper's speed-ratio boundary): in its canonical geometry # r=1, pursuer=(0,y), the tangent arc has beta=2 atan(y/r), while the # pursuer tangent distance is y. Evader reaches S first iff # rho=ve/vp >= beta*r/y; this threshold is <=2 and tends to 2 as y->0. ys = np.geomspace(1e-4, 20.0, 80) ratio_thresholds = (2*np.arctan(ys)/ys) speed_sweep = { "y_over_r": ys.tolist(), "observed_threshold_ve_over_vp": ratio_thresholds.tolist(), "predicted_supremum": 2.0, "observed_supremum": float(np.max(ratio_thresholds)), "threshold_at_y_over_r_1": float(ratio_thresholds[np.argmin(abs(ys-1.0))]), } # Tiny comparison: random residual field, projected layer versus raw policy. def raw_fn(x): return np.array([-2., .35*np.sin(3*x[1])]) def safe_fn(x): return project_velocity(x, raw_fn(x), c, r, lookahead_dt=.005) raw_end, raw_min=integrate(x0,raw_fn,c,r) safe_end, safe_min=integrate(x0,safe_fn,c,r) result={ "math_check": {"max_relative_tangent_distance_error":float(max(rel_errors)), "max_tangent_clearance_error":float(max(clear_errors))}, "branch_switch": {"predicted_heading_rad":0.0,"observed_heading_rad":float(switch), "max_time_difference_at_switch":float(abs(diffs[np.argmin(abs(diffs))]))}, "safety_sweep": {"residual_magnitudes":np.linspace(.1,8.,12).tolist(), "raw_min_clearance":raw_clear,"projected_min_clearance":projected_clear}, "speed_ratio_sweep": speed_sweep, "mini_experiment": {"raw_min_clearance":raw_min,"projected_min_clearance":safe_min, "raw_final_position":raw_end.tolist(),"projected_final_position":safe_end.tolist()}, } Path('results.json').write_text(json.dumps(result, indent=2)) print(json.dumps(result, indent=2)) if __name__ == '__main__': main()