1import json, math
 2import numpy as np
 3from scipy.linalg import solve_discrete_lyapunov, eigvalsh
 4from event_triggered_preconditioner import make_hessians, lmi_margin
 5
 6
 7def fixed_run(H, K, h, sigma, steps=80, max_interval=10000):
 8    x=np.array([2.0,-1.5]); last=x.copy(); events=0; since=0; norms=[]
 9    for t in range(steps):
10        x=x-h*K@(H@last); since+=1
11        d=last-x
12        if d@d > sigma*sigma*(x@x+1e-12) or since>=max_interval:
13            last=x.copy(); events+=1; since=0
14        norms.append(np.linalg.norm(x))
15        if not np.isfinite(norms[-1]) or norms[-1]>1e10: return False, events, norms
16    return True,events,norms
17
18
19def main():
20    Hs=make_hessians(); H=Hs[0]
21    Hbar=sum(Hs)/2
22    Ks={'diag':np.diag(.72/np.diag(Hbar)), 'dense':.72*np.linalg.inv(Hbar)}
23    print('PREDICTION_1 discrete fixed-H boundary: h_boundary=2/rho(KH)')
24    for name,K in Ks.items():
25        rho=max(abs(np.linalg.eigvals(K@H)))
26        hb=2/rho
27        rows=[]
28        for mult in [.90,.99,1.01,1.10]:
29            ok,_,norms=fixed_run(H,K,mult*hb,0,steps=100,max_interval=1)
30            rows.append({'mult':mult,'h':mult*hb,'stable_100_steps':ok,'final_norm':float(norms[-1])})
31        print(json.dumps({'method':name,'rho':float(rho),'predicted_boundary':float(hb),'observed':rows}))
32
33    print('PREDICTION_2 LMI threshold: margin crosses zero as sigma increases')
34    for name,K in Ks.items():
35        lo,hi=0.,10.
36        # Fixed P and lambda as in the implementation; locate the numerical certificate boundary.
37        while lmi_margin(Hs,K,.02,hi,lam=1)[0]<0: hi*=2
38        for _ in range(60):
39            mid=(lo+hi)/2
40            if lmi_margin(Hs,K,.02,mid,lam=1)[0]<0: lo=mid
41            else: hi=mid
42        vals=[]
43        for s in [0,.25*lo,.5*lo,.75*lo,lo,1.05*lo]:
44            mar,_=lmi_margin(Hs,K,.02,s,lam=1)
45            vals.append({'sigma':float(s),'margin':float(mar),'certified':bool(mar<0)})
46        print(json.dumps({'method':name,'predicted_sigma_boundary':float(lo),'observed':vals}))
47
48    print('PREDICTION_3 finite-horizon trigger scaling: larger sigma lowers events, with bounded error')
49    H=Hs[0]; K=Ks['dense']; h=.25
50    # This h is safely below the fixed-H synchronous boundary and exposes stale updates.
51    for sigma in [0,.02,.05,.10,.20,.40]:
52        ok,ev,norms=fixed_run(H,K,h,sigma,steps=40,max_interval=1000)
53        sync=fixed_run(H,K,h,0,steps=40,max_interval=1)[2][-1]
54        print(json.dumps({'sigma':sigma,'events':ev,'event_fraction':ev/40,
55            'norm_at_40':float(norms[-1]),'ratio_to_sync_norm':float(norms[-1]/(sync+1e-30)),
56            'stable':ok}))
57
58if __name__=='__main__': main()